From d478435b01bca4a858ddb194e5d172f47292f804 Mon Sep 17 00:00:00 2001 From: "copilot-swe-agent[bot]" <198982749+Copilot@users.noreply.github.com> Date: Thu, 22 May 2025 20:21:33 +0000 Subject: [PATCH 01/21] Initial plan for issue From 931fcaaf7342a668810771deb71d25483cc8afae Mon Sep 17 00:00:00 2001 From: "copilot-swe-agent[bot]" <198982749+Copilot@users.noreply.github.com> Date: Thu, 22 May 2025 20:24:19 +0000 Subject: [PATCH 02/21] Create template matching reproducibility test script Co-authored-by: bHimes <3077528+bHimes@users.noreply.github.com> --- .../test_template_reproducibility.py | 96 +++++++++++++++++++ 1 file changed, 96 insertions(+) create mode 100755 scripts/testing/programs/match_template/test_template_reproducibility.py diff --git a/scripts/testing/programs/match_template/test_template_reproducibility.py b/scripts/testing/programs/match_template/test_template_reproducibility.py new file mode 100755 index 000000000..04373fa51 --- /dev/null +++ b/scripts/testing/programs/match_template/test_template_reproducibility.py @@ -0,0 +1,96 @@ +#!/usr/bin/env python3 +import annoying_hack + +from os.path import join +from os import makedirs +import cistem_test_utils.args as tmArgs +import cistem_test_utils.make_tmp_runfile as mktmp +import cistem_test_utils.run_job as runner +import mrcfile +import numpy as np +import tempfile +import os + +# By default the "_gpu" suffix will be added unless the --old-cistem flag is used +# or the --cpu flag is used +wanted_binary_name = 'match_template' + + +def main(): + args = tmArgs.parse_TM_args(wanted_binary_name) + + # Create a temporary directory to store our replicate MIPs + temp_dir = tempfile.mkdtemp(dir="/tmp", prefix="template_match_reproducibility_") + print(f"Temporary directory created at: {temp_dir}") + + # We'll run template matching 3 times to generate replicates + elapsed_time = [0, 0, 0] + mip_filenames = [] + + # Run the template matching 3 times + for replicate in range(0, 3): + # Use Apoferritin dataset with image 0 + config = tmArgs.get_config(args, 'Apoferritin', 0, 0) + + # Run the template matching + tmp_filename_match_template, tmp_filename_make_template_results = mktmp.make_tmp_runfile(config) + elapsed_time[replicate] = runner.run_job(tmp_filename_match_template) + runner.run_job(tmp_filename_make_template_results) + + # Copy the MIP file to our temp directory with a unique name + mip_src = join(config['output_file_prefix'], 'mip.mrc') + mip_dst = join(temp_dir, f"mip_replicate_{replicate+1}.mrc") + + # Use os.system to copy the file - using cp preserves the MRC file format + os.system(f"cp {mip_src} {mip_dst}") + mip_filenames.append(mip_dst) + + print(f"Completed replicate {replicate+1}, time: {elapsed_time[replicate]:.2f}s") + + # Now load the MIP files and compare them + mip_data = [] + for filename in mip_filenames: + with mrcfile.open(filename) as mrc: + mip_data.append(mrc.data) + + # Compute similarity metrics between replicates + print("\nReproducibility Analysis:") + print("========================") + + # Pairwise comparisons + pairs = [(0,1), (0,2), (1,2)] + for i, j in pairs: + # Calculate correlation coefficient + correlation = np.corrcoef(mip_data[i].flatten(), mip_data[j].flatten())[0,1] + + # Calculate mean absolute difference + mean_abs_diff = np.mean(np.abs(mip_data[i] - mip_data[j])) + + # Calculate peak signal to noise ratio (PSNR) + mse = np.mean((mip_data[i] - mip_data[j]) ** 2) + max_pixel = max(np.max(mip_data[i]), np.max(mip_data[j])) + psnr = 20 * np.log10(max_pixel) - 10 * np.log10(mse) if mse > 0 else float('inf') + + print(f"Comparing replicate {i+1} vs {j+1}:") + print(f" Correlation coefficient: {correlation:.6f}") + print(f" Mean absolute difference: {mean_abs_diff:.6f}") + print(f" Peak signal-to-noise ratio (PSNR): {psnr:.2f} dB") + + # Overall similarity across all three replicates + all_correlations = [] + for i, j in pairs: + corr = np.corrcoef(mip_data[i].flatten(), mip_data[j].flatten())[0,1] + all_correlations.append(corr) + + print("\nOverall reproducibility:") + print(f" Mean correlation: {np.mean(all_correlations):.6f}") + print(f" Min correlation: {np.min(all_correlations):.6f}") + print(f" Max correlation: {np.max(all_correlations):.6f}") + + # Print the directory where files are saved + print(f"\nMIP files saved in: {temp_dir}") + + +# Check if main function and run +if __name__ == '__main__': + main() \ No newline at end of file From a8ca4068ac7e369163dadabc26e3d74be3f9b27c Mon Sep 17 00:00:00 2001 From: "copilot-swe-agent[bot]" <198982749+Copilot@users.noreply.github.com> Date: Thu, 22 May 2025 20:25:16 +0000 Subject: [PATCH 03/21] Improve error handling and documentation in template matching test Co-authored-by: bHimes <3077528+bHimes@users.noreply.github.com> --- .../test_template_reproducibility.py | 198 ++++++++++++------ 1 file changed, 132 insertions(+), 66 deletions(-) diff --git a/scripts/testing/programs/match_template/test_template_reproducibility.py b/scripts/testing/programs/match_template/test_template_reproducibility.py index 04373fa51..2637c3aae 100755 --- a/scripts/testing/programs/match_template/test_template_reproducibility.py +++ b/scripts/testing/programs/match_template/test_template_reproducibility.py @@ -1,7 +1,20 @@ #!/usr/bin/env python3 -import annoying_hack +""" +Template Matching Reproducibility Test + +This script runs the cisTEM template matching GPU binary three times using the same input data +and parameters, then analyzes the reproducibility by comparing the resulting MIP images. + +It uses the Apoferritin dataset for testing, and saves the MIP (Maximum Intensity Projection) +files to a temporary directory with unique filenames for each replicate. -from os.path import join +The script loads the MIP images using the mrcfile package and compares them using numpy to +measure pixel similarity between replicates with various metrics. +""" + +import annoying_hack +import sys +from os.path import join, exists from os import makedirs import cistem_test_utils.args as tmArgs import cistem_test_utils.make_tmp_runfile as mktmp @@ -10,6 +23,8 @@ import numpy as np import tempfile import os +import shutil +import subprocess # By default the "_gpu" suffix will be added unless the --old-cistem flag is used # or the --cpu flag is used @@ -17,80 +32,131 @@ def main(): - args = tmArgs.parse_TM_args(wanted_binary_name) + try: + # Parse command-line arguments + args = tmArgs.parse_TM_args(wanted_binary_name) - # Create a temporary directory to store our replicate MIPs - temp_dir = tempfile.mkdtemp(dir="/tmp", prefix="template_match_reproducibility_") - print(f"Temporary directory created at: {temp_dir}") - - # We'll run template matching 3 times to generate replicates - elapsed_time = [0, 0, 0] - mip_filenames = [] - - # Run the template matching 3 times - for replicate in range(0, 3): - # Use Apoferritin dataset with image 0 - config = tmArgs.get_config(args, 'Apoferritin', 0, 0) + # Create a temporary directory to store our replicate MIPs + temp_dir = tempfile.mkdtemp(dir="/tmp", prefix="template_match_reproducibility_") + print(f"Temporary directory created at: {temp_dir}") - # Run the template matching - tmp_filename_match_template, tmp_filename_make_template_results = mktmp.make_tmp_runfile(config) - elapsed_time[replicate] = runner.run_job(tmp_filename_match_template) - runner.run_job(tmp_filename_make_template_results) + # We'll run template matching 3 times to generate replicates + elapsed_time = [0, 0, 0] + mip_filenames = [] - # Copy the MIP file to our temp directory with a unique name - mip_src = join(config['output_file_prefix'], 'mip.mrc') - mip_dst = join(temp_dir, f"mip_replicate_{replicate+1}.mrc") + # Run the template matching 3 times + for replicate in range(0, 3): + try: + # Use Apoferritin dataset with image 0 + config = tmArgs.get_config(args, 'Apoferritin', 0, 0) + + # Run the template matching + tmp_filename_match_template, tmp_filename_make_template_results = mktmp.make_tmp_runfile(config) + elapsed_time[replicate] = runner.run_job(tmp_filename_match_template) + runner.run_job(tmp_filename_make_template_results) + + # Copy the MIP file to our temp directory with a unique name + mip_src = join(config['output_file_prefix'], 'mip.mrc') + + # Check if source file exists + if not exists(mip_src): + raise FileNotFoundError(f"MIP file not found at {mip_src}") + + mip_dst = join(temp_dir, f"mip_replicate_{replicate+1}.mrc") + + # Use copy2 from shutil to preserve metadata + shutil.copy2(mip_src, mip_dst) + mip_filenames.append(mip_dst) + + print(f"Completed replicate {replicate+1}, time: {elapsed_time[replicate]:.2f}s") + + except Exception as e: + print(f"Error during replicate {replicate+1}: {str(e)}") + continue - # Use os.system to copy the file - using cp preserves the MRC file format - os.system(f"cp {mip_src} {mip_dst}") - mip_filenames.append(mip_dst) + # Check if we have all three replicates + if len(mip_filenames) != 3: + print(f"Warning: Only {len(mip_filenames)} replicates were successfully processed") + if len(mip_filenames) < 2: + print("At least 2 replicates are required for comparison. Exiting.") + return 1 - print(f"Completed replicate {replicate+1}, time: {elapsed_time[replicate]:.2f}s") - - # Now load the MIP files and compare them - mip_data = [] - for filename in mip_filenames: - with mrcfile.open(filename) as mrc: - mip_data.append(mrc.data) - - # Compute similarity metrics between replicates - print("\nReproducibility Analysis:") - print("========================") - - # Pairwise comparisons - pairs = [(0,1), (0,2), (1,2)] - for i, j in pairs: - # Calculate correlation coefficient - correlation = np.corrcoef(mip_data[i].flatten(), mip_data[j].flatten())[0,1] + # Now load the MIP files and compare them + print("\nLoading MIP files for analysis...") + mip_data = [] + for filename in mip_filenames: + try: + with mrcfile.open(filename) as mrc: + mip_data.append(mrc.data) + print(f"Loaded {filename}, shape: {mrc.data.shape}, dtype: {mrc.data.dtype}") + except Exception as e: + print(f"Error loading {filename}: {str(e)}") - # Calculate mean absolute difference - mean_abs_diff = np.mean(np.abs(mip_data[i] - mip_data[j])) + # Check if we have at least 2 MIP files loaded + if len(mip_data) < 2: + print("Failed to load at least 2 MIP files. Cannot perform comparison.") + return 1 - # Calculate peak signal to noise ratio (PSNR) - mse = np.mean((mip_data[i] - mip_data[j]) ** 2) - max_pixel = max(np.max(mip_data[i]), np.max(mip_data[j])) - psnr = 20 * np.log10(max_pixel) - 10 * np.log10(mse) if mse > 0 else float('inf') + # Compute similarity metrics between replicates + print("\nReproducibility Analysis:") + print("========================") - print(f"Comparing replicate {i+1} vs {j+1}:") - print(f" Correlation coefficient: {correlation:.6f}") - print(f" Mean absolute difference: {mean_abs_diff:.6f}") - print(f" Peak signal-to-noise ratio (PSNR): {psnr:.2f} dB") - - # Overall similarity across all three replicates - all_correlations = [] - for i, j in pairs: - corr = np.corrcoef(mip_data[i].flatten(), mip_data[j].flatten())[0,1] - all_correlations.append(corr) - - print("\nOverall reproducibility:") - print(f" Mean correlation: {np.mean(all_correlations):.6f}") - print(f" Min correlation: {np.min(all_correlations):.6f}") - print(f" Max correlation: {np.max(all_correlations):.6f}") + # Pairwise comparisons of available replicates + pairs = [(i, j) for i in range(len(mip_data)) for j in range(i+1, len(mip_data))] + + all_correlations = [] + all_mean_abs_diffs = [] + all_psnrs = [] + + for i, j in pairs: + try: + # Calculate correlation coefficient + correlation = np.corrcoef(mip_data[i].flatten(), mip_data[j].flatten())[0,1] + all_correlations.append(correlation) + + # Calculate mean absolute difference + mean_abs_diff = np.mean(np.abs(mip_data[i] - mip_data[j])) + all_mean_abs_diffs.append(mean_abs_diff) + + # Calculate peak signal to noise ratio (PSNR) + mse = np.mean((mip_data[i] - mip_data[j]) ** 2) + max_pixel = max(np.max(mip_data[i]), np.max(mip_data[j])) + psnr = 20 * np.log10(max_pixel) - 10 * np.log10(mse) if mse > 0 else float('inf') + all_psnrs.append(psnr) + + # Print results for this pair + print(f"Comparing replicate {i+1} vs {j+1}:") + print(f" Correlation coefficient: {correlation:.6f}") + print(f" Mean absolute difference: {mean_abs_diff:.6f}") + print(f" Peak signal-to-noise ratio (PSNR): {psnr:.2f} dB") + + except Exception as e: + print(f"Error comparing replicates {i+1} and {j+1}: {str(e)}") + + # Overall similarity across all replicates + if all_correlations: + print("\nOverall reproducibility:") + print(f" Mean correlation: {np.mean(all_correlations):.6f}") + print(f" Min correlation: {np.min(all_correlations):.6f}") + print(f" Max correlation: {np.max(all_correlations):.6f}") + print(f" Mean absolute diff (avg): {np.mean(all_mean_abs_diffs):.6f}") + print(f" PSNR (avg): {np.mean(all_psnrs):.2f} dB") + + # Print the directory where files are saved + print(f"\nMIP files saved in: {temp_dir}") + + return 0 - # Print the directory where files are saved - print(f"\nMIP files saved in: {temp_dir}") + except KeyboardInterrupt: + print("\nTest interrupted by user.") + return 130 + except Exception as e: + print(f"Error in template matching reproducibility test: {str(e)}") + import traceback + traceback.print_exc() + return 1 # Check if main function and run if __name__ == '__main__': - main() \ No newline at end of file + sys.exit(main()) \ No newline at end of file From 1a3b45637f89ce747ddf4f71477a622a35085843 Mon Sep 17 00:00:00 2001 From: "copilot-swe-agent[bot]" <198982749+Copilot@users.noreply.github.com> Date: Thu, 22 May 2025 20:27:12 +0000 Subject: [PATCH 04/21] Optimize template matching reproducibility test Co-authored-by: bHimes <3077528+bHimes@users.noreply.github.com> --- .../test_template_reproducibility.py | 26 +++++++++---------- 1 file changed, 13 insertions(+), 13 deletions(-) diff --git a/scripts/testing/programs/match_template/test_template_reproducibility.py b/scripts/testing/programs/match_template/test_template_reproducibility.py index 2637c3aae..9e981b752 100755 --- a/scripts/testing/programs/match_template/test_template_reproducibility.py +++ b/scripts/testing/programs/match_template/test_template_reproducibility.py @@ -13,7 +13,7 @@ """ import annoying_hack -import sys + from os.path import join, exists from os import makedirs import cistem_test_utils.args as tmArgs @@ -24,7 +24,7 @@ import tempfile import os import shutil -import subprocess +import sys # By default the "_gpu" suffix will be added unless the --old-cistem flag is used # or the --cpu flag is used @@ -50,24 +50,24 @@ def main(): # Use Apoferritin dataset with image 0 config = tmArgs.get_config(args, 'Apoferritin', 0, 0) + # Create a unique output file prefix for each replicate + original_prefix = config['output_file_prefix'] + config['output_file_prefix'] = join(temp_dir, f"replicate_{replicate+1}") + makedirs(config['output_file_prefix'], exist_ok=True) + # Run the template matching tmp_filename_match_template, tmp_filename_make_template_results = mktmp.make_tmp_runfile(config) elapsed_time[replicate] = runner.run_job(tmp_filename_match_template) runner.run_job(tmp_filename_make_template_results) - # Copy the MIP file to our temp directory with a unique name - mip_src = join(config['output_file_prefix'], 'mip.mrc') - - # Check if source file exists - if not exists(mip_src): - raise FileNotFoundError(f"MIP file not found at {mip_src}") - - mip_dst = join(temp_dir, f"mip_replicate_{replicate+1}.mrc") + # The MIP file is already in our temp directory with the unique replicate prefix + mip_file = join(config['output_file_prefix'], 'mip.mrc') - # Use copy2 from shutil to preserve metadata - shutil.copy2(mip_src, mip_dst) - mip_filenames.append(mip_dst) + # Check if MIP file exists + if not exists(mip_file): + raise FileNotFoundError(f"MIP file not found at {mip_file}") + mip_filenames.append(mip_file) print(f"Completed replicate {replicate+1}, time: {elapsed_time[replicate]:.2f}s") except Exception as e: From 6548b4a012628cb2221589eb2960857b6d4e7633 Mon Sep 17 00:00:00 2001 From: "copilot-swe-agent[bot]" <198982749+Copilot@users.noreply.github.com> Date: Thu, 22 May 2025 20:51:30 +0000 Subject: [PATCH 05/21] Add --fast_fft flag to template matching tests Co-authored-by: bHimes <3077528+bHimes@users.noreply.github.com> --- scripts/testing/programs/cistem_test_utils/args.py | 4 ++++ .../programs/cistem_test_utils/make_tmp_runfile.py | 12 ++++++++++++ 2 files changed, 16 insertions(+) diff --git a/scripts/testing/programs/cistem_test_utils/args.py b/scripts/testing/programs/cistem_test_utils/args.py index 573c2c8b1..4a0ec6c01 100644 --- a/scripts/testing/programs/cistem_test_utils/args.py +++ b/scripts/testing/programs/cistem_test_utils/args.py @@ -107,6 +107,10 @@ def parse_TM_args(wanted_binary_name): help='Use this flag if you are using the cpu version of cisTEM') args_to_check.append('cpu') + parser.add_argument('--fast_fft', action='store_true', default=True, + help='Use fast FFT implementation (default: True)') + args_to_check.append('fast_fft') + args = parser.parse_args() args.binary_name = wanted_binary_name diff --git a/scripts/testing/programs/cistem_test_utils/make_tmp_runfile.py b/scripts/testing/programs/cistem_test_utils/make_tmp_runfile.py index 687d07869..0e61e1802 100644 --- a/scripts/testing/programs/cistem_test_utils/make_tmp_runfile.py +++ b/scripts/testing/programs/cistem_test_utils/make_tmp_runfile.py @@ -9,6 +9,17 @@ def match_template(config): use_gpu = "no" else: use_gpu = "yes" + + # Check if fast_fft is requested but GPU is not enabled + if config.get('fast_fft') and use_gpu == "no": + print("Error: --fast_fft option requires GPU to be enabled. Cannot use fast FFT with CPU.") + exit(1) + + # Set the use_fast_fft parameter based on the fast_fft flag + if config.get('fast_fft'): + use_fast_fft = "yes" + else: + use_fast_fft = "no" high_res_limit = 2*config.get('data')[config.get('img_number')].get('pixel_size') @@ -46,6 +57,7 @@ def match_template(config): str(config.get('mask_radius')), str(config.get('model')[config.get('ref_number')].get('symmetry')), use_gpu, + use_fast_fft, str(config.get('max_threads'))] pre_process_cmd = " " From 9ac268081cc353d4e76d0cff2a6752379c2f9495 Mon Sep 17 00:00:00 2001 From: "copilot-swe-agent[bot]" <198982749+Copilot@users.noreply.github.com> Date: Thu, 22 May 2025 21:11:59 +0000 Subject: [PATCH 06/21] Implement temp directory management and update FastFFT help text Co-authored-by: bHimes <3077528+bHimes@users.noreply.github.com> --- .cistem_temp_dirs.log | 7 + .../programs/cistem_test_utils/__init__.py | 6 + .../programs/cistem_test_utils/args.py | 12 +- .../cistem_test_utils/temp_dir_manager.py | 172 ++++++++++++++++++ .../test_template_reproducibility.py | 25 ++- 5 files changed, 219 insertions(+), 3 deletions(-) create mode 100644 .cistem_temp_dirs.log create mode 100644 scripts/testing/programs/cistem_test_utils/temp_dir_manager.py diff --git a/.cistem_temp_dirs.log b/.cistem_temp_dirs.log new file mode 100644 index 000000000..3c6d38ec9 --- /dev/null +++ b/.cistem_temp_dirs.log @@ -0,0 +1,7 @@ +[ + { + "path": "/tmp/cistem_test_ndv2z151", + "created_at": "2025-05-22T21:11:37.778761", + "key": "cistem_test_ndv2z151" + } +] \ No newline at end of file diff --git a/scripts/testing/programs/cistem_test_utils/__init__.py b/scripts/testing/programs/cistem_test_utils/__init__.py index e69de29bb..756dbcee7 100644 --- a/scripts/testing/programs/cistem_test_utils/__init__.py +++ b/scripts/testing/programs/cistem_test_utils/__init__.py @@ -0,0 +1,6 @@ +# cisTEM testing utilities package +# Import commonly used modules to make them easily accessible +from . import args +from . import make_tmp_runfile +from . import run_job +from . import temp_dir_manager diff --git a/scripts/testing/programs/cistem_test_utils/args.py b/scripts/testing/programs/cistem_test_utils/args.py index 4a0ec6c01..cad2fc731 100644 --- a/scripts/testing/programs/cistem_test_utils/args.py +++ b/scripts/testing/programs/cistem_test_utils/args.py @@ -77,6 +77,16 @@ def parse_TM_args(wanted_binary_name): # To help check, we'll keep this dict to use as a check. args_to_check = [] parser = argparse.ArgumentParser(description='Test the k3 rotation binary') + + # Temp directory management arguments + temp_management_group = parser.add_argument_group('Temporary Directory Management') + temp_management_group.add_argument('--list-temp-dirs', action='store_true', + help='List all tracked temporary directories') + temp_management_group.add_argument('--rm-temp-dir', type=int, metavar='INDEX', + help='Remove a specific temporary directory by index') + temp_management_group.add_argument('--rm-all-temp-dirs', action='store_true', + help='Remove all tracked temporary directories') + # Required argument parser.add_argument( '--binary_path', help='Path to the directory with the binary to be tested (Required)', required=True) @@ -108,7 +118,7 @@ def parse_TM_args(wanted_binary_name): args_to_check.append('cpu') parser.add_argument('--fast_fft', action='store_true', default=True, - help='Use fast FFT implementation (default: True)') + help='Use FastFFT implementation (default: True)') args_to_check.append('fast_fft') args = parser.parse_args() diff --git a/scripts/testing/programs/cistem_test_utils/temp_dir_manager.py b/scripts/testing/programs/cistem_test_utils/temp_dir_manager.py new file mode 100644 index 000000000..eb6c69c5d --- /dev/null +++ b/scripts/testing/programs/cistem_test_utils/temp_dir_manager.py @@ -0,0 +1,172 @@ +""" +Temporary Directory Manager for cisTEM tests. + +This module provides functionality to track, list, and clean up temporary directories +created during cisTEM test runs. It maintains a log file in the current working directory +that records information about each temporary directory including its path, creation time, +and a unique identifier. +""" + +import os +import tempfile +import datetime +import json +import shutil +from typing import List, Dict, Optional, Tuple, Union + +# The log file will be a hidden file in the current working directory +LOG_FILE = os.path.join(os.getcwd(), '.cistem_temp_dirs.log') + +def create_temp_dir(prefix: str = "cistem_test_", dir: str = "/tmp") -> str: + """ + Create a temporary directory and log its information. + + Args: + prefix: Prefix for the temporary directory name + dir: Parent directory where the temp directory will be created + + Returns: + Path to the created temporary directory + """ + # Create the temporary directory + temp_dir = tempfile.mkdtemp(dir=dir, prefix=prefix) + + # Log the created directory + log_entry = { + 'path': temp_dir, + 'created_at': datetime.datetime.now().isoformat(), + 'key': os.path.basename(temp_dir) # Use the basename as the unique key + } + + # Read existing log if it exists + log_entries = [] + if os.path.exists(LOG_FILE): + try: + with open(LOG_FILE, 'r') as f: + log_entries = json.load(f) + except json.JSONDecodeError: + # If the file is corrupted, start with an empty list + log_entries = [] + + # Add the new entry and write back to the log file + log_entries.append(log_entry) + with open(LOG_FILE, 'w') as f: + json.dump(log_entries, f, indent=2) + + return temp_dir + +def list_temp_dirs() -> List[Dict]: + """ + List all logged temporary directories. + + Returns: + List of dictionaries containing information about each temp directory + """ + if not os.path.exists(LOG_FILE): + return [] + + try: + with open(LOG_FILE, 'r') as f: + log_entries = json.load(f) + + # Filter out directories that no longer exist + valid_entries = [] + for entry in log_entries: + if os.path.exists(entry['path']): + valid_entries.append(entry) + + # Update the log file if some directories were already removed + if len(valid_entries) != len(log_entries): + with open(LOG_FILE, 'w') as f: + json.dump(valid_entries, f, indent=2) + + return valid_entries + except json.JSONDecodeError: + # If the file is corrupted, return an empty list + return [] + +def remove_temp_dir(index: int) -> Tuple[bool, str]: + """ + Remove a temporary directory by its index in the log. + + Args: + index: The index of the directory to remove (0-based) + + Returns: + Tuple of (success, message) + """ + log_entries = list_temp_dirs() + + if not log_entries: + return False, "No temporary directories found in the log." + + if index < 0 or index >= len(log_entries): + return False, f"Invalid index: {index}. Valid range is 0-{len(log_entries)-1}." + + entry = log_entries[index] + path = entry['path'] + + # Check if the directory exists + if not os.path.exists(path): + # Remove from log and return + log_entries.pop(index) + with open(LOG_FILE, 'w') as f: + json.dump(log_entries, f, indent=2) + return False, f"Directory {path} no longer exists. Removed from log." + + # Remove the directory + try: + shutil.rmtree(path) + # Update the log + log_entries.pop(index) + with open(LOG_FILE, 'w') as f: + json.dump(log_entries, f, indent=2) + return True, f"Successfully removed {path}" + except Exception as e: + return False, f"Failed to remove {path}: {str(e)}" + +def remove_all_temp_dirs() -> Tuple[int, int]: + """ + Remove all temporary directories in the log. + + Returns: + Tuple of (success_count, failure_count) + """ + log_entries = list_temp_dirs() + success_count = 0 + failure_count = 0 + + for entry in log_entries: + path = entry['path'] + if os.path.exists(path): + try: + shutil.rmtree(path) + success_count += 1 + except: + failure_count += 1 + + # Clear the log file + with open(LOG_FILE, 'w') as f: + json.dump([], f) + + return success_count, failure_count + +def print_temp_dirs() -> None: + """ + Print information about all logged temporary directories. + """ + log_entries = list_temp_dirs() + + if not log_entries: + print("No temporary directories found.") + return + + print(f"Found {len(log_entries)} temporary directories:") + print("-" * 80) + for i, entry in enumerate(log_entries): + created_at = datetime.datetime.fromisoformat(entry['created_at']) + formatted_time = created_at.strftime("%Y-%m-%d %H:%M:%S") + print(f"[{i}] {entry['path']}") + print(f" Created: {formatted_time}") + print(f" Key: {entry['key']}") + print("-" * 80) \ No newline at end of file diff --git a/scripts/testing/programs/match_template/test_template_reproducibility.py b/scripts/testing/programs/match_template/test_template_reproducibility.py index 9e981b752..3349877ca 100755 --- a/scripts/testing/programs/match_template/test_template_reproducibility.py +++ b/scripts/testing/programs/match_template/test_template_reproducibility.py @@ -19,6 +19,7 @@ import cistem_test_utils.args as tmArgs import cistem_test_utils.make_tmp_runfile as mktmp import cistem_test_utils.run_job as runner +import cistem_test_utils.temp_dir_manager as temp_dir_manager import mrcfile import numpy as np import tempfile @@ -35,9 +36,26 @@ def main(): try: # Parse command-line arguments args = tmArgs.parse_TM_args(wanted_binary_name) + + # Handle temp directory management options + if args.list_temp_dirs: + temp_dir_manager.print_temp_dirs() + return 0 + + if args.rm_temp_dir is not None: + success, message = temp_dir_manager.remove_temp_dir(args.rm_temp_dir) + print(message) + return 0 if success else 1 + + if args.rm_all_temp_dirs: + success_count, failure_count = temp_dir_manager.remove_all_temp_dirs() + print(f"Successfully removed {success_count} temporary directories.") + if failure_count > 0: + print(f"Failed to remove {failure_count} temporary directories.") + return 0 if failure_count == 0 else 1 - # Create a temporary directory to store our replicate MIPs - temp_dir = tempfile.mkdtemp(dir="/tmp", prefix="template_match_reproducibility_") + # Create a temporary directory to store our replicate MIPs and track it + temp_dir = temp_dir_manager.create_temp_dir(prefix="template_match_reproducibility_") print(f"Temporary directory created at: {temp_dir}") # We'll run template matching 3 times to generate replicates @@ -144,6 +162,9 @@ def main(): # Print the directory where files are saved print(f"\nMIP files saved in: {temp_dir}") + print("To list temp directories: --list-temp-dirs") + print("To remove this directory: --rm-temp-dir INDEX") + print("To remove all temp directories: --rm-all-temp-dirs") return 0 From b3c7f6028dc368e1bea5105f4279832ab4fc62b4 Mon Sep 17 00:00:00 2001 From: "copilot-swe-agent[bot]" <198982749+Copilot@users.noreply.github.com> Date: Thu, 22 May 2025 21:19:08 +0000 Subject: [PATCH 07/21] Add threshold extraction and relative error reporting to template matching test Co-authored-by: bHimes <3077528+bHimes@users.noreply.github.com> --- .../test_template_reproducibility.py | 80 +++++++++++++++++++ 1 file changed, 80 insertions(+) diff --git a/scripts/testing/programs/match_template/test_template_reproducibility.py b/scripts/testing/programs/match_template/test_template_reproducibility.py index 3349877ca..dad7a8630 100755 --- a/scripts/testing/programs/match_template/test_template_reproducibility.py +++ b/scripts/testing/programs/match_template/test_template_reproducibility.py @@ -26,12 +26,46 @@ import os import shutil import sys +import re # By default the "_gpu" suffix will be added unless the --old-cistem flag is used # or the --cpu flag is used wanted_binary_name = 'match_template' +def extract_threshold_value(hist_file_path): + """ + Extract the threshold value from the histogram text file. + + The threshold value is in the first line of the file, which starts with + "# Expected". This function extracts the numerical value following this marker. + + Args: + hist_file_path (str): Path to the histogram text file + + Returns: + float: The extracted threshold value + + Raises: + FileNotFoundError: If the file doesn't exist + ValueError: If the threshold value cannot be found or parsed + """ + if not exists(hist_file_path): + raise FileNotFoundError(f"Histogram file not found at {hist_file_path}") + + with open(hist_file_path, 'r') as f: + first_line = f.readline().strip() + + # Use regex to extract the threshold value from the line + # Looking for a pattern like "# Expected max is: 0.12345" + match = re.search(r'# Expected max is: ([\d.e+-]+)', first_line) + if not match: + raise ValueError(f"Could not find threshold value in {hist_file_path}") + + threshold_value = float(match.group(1)) + return threshold_value + + def main(): try: # Parse command-line arguments @@ -61,6 +95,8 @@ def main(): # We'll run template matching 3 times to generate replicates elapsed_time = [0, 0, 0] mip_filenames = [] + hist_filenames = [] + threshold_values = [] # Run the template matching 3 times for replicate in range(0, 3): @@ -85,8 +121,16 @@ def main(): if not exists(mip_file): raise FileNotFoundError(f"MIP file not found at {mip_file}") + # Get the histogram file path and extract the threshold value + hist_file = join(config['output_file_prefix'], 'hist.txt') + threshold_value = extract_threshold_value(hist_file) + mip_filenames.append(mip_file) + hist_filenames.append(hist_file) + threshold_values.append(threshold_value) + print(f"Completed replicate {replicate+1}, time: {elapsed_time[replicate]:.2f}s") + print(f" Threshold value: {threshold_value:.6e}") except Exception as e: print(f"Error during replicate {replicate+1}: {str(e)}") @@ -114,6 +158,20 @@ def main(): if len(mip_data) < 2: print("Failed to load at least 2 MIP files. Cannot perform comparison.") return 1 + + # Verify all threshold values are the same + if threshold_values: + if not all(abs(v - threshold_values[0]) < 1e-10 for v in threshold_values): + print("Warning: Threshold values differ between replicates:") + for i, val in enumerate(threshold_values): + print(f" Replicate {i+1}: {val:.6e}") + + # Use the first threshold value for calculations + threshold_value = threshold_values[0] + print(f"\nThreshold value: {threshold_value:.6e}") + else: + print("Warning: No threshold values were extracted.") + threshold_value = None # Compute similarity metrics between replicates print("\nReproducibility Analysis:") @@ -136,6 +194,14 @@ def main(): mean_abs_diff = np.mean(np.abs(mip_data[i] - mip_data[j])) all_mean_abs_diffs.append(mean_abs_diff) + # Calculate relative error if threshold value is available + if threshold_value and threshold_value > 0: + relative_error = mean_abs_diff / threshold_value + relative_error_ppm = relative_error * 1e6 # Parts per million + else: + relative_error = None + relative_error_ppm = None + # Calculate peak signal to noise ratio (PSNR) mse = np.mean((mip_data[i] - mip_data[j]) ** 2) max_pixel = max(np.max(mip_data[i]), np.max(mip_data[j])) @@ -146,6 +212,12 @@ def main(): print(f"Comparing replicate {i+1} vs {j+1}:") print(f" Correlation coefficient: {correlation:.6f}") print(f" Mean absolute difference: {mean_abs_diff:.6f}") + + # Print relative error if available + if relative_error is not None: + print(f" Relative error (scientific): {relative_error:.6e}") + print(f" Relative error (ppm): {relative_error_ppm:.2f} ppm") + print(f" Peak signal-to-noise ratio (PSNR): {psnr:.2f} dB") except Exception as e: @@ -158,6 +230,14 @@ def main(): print(f" Min correlation: {np.min(all_correlations):.6f}") print(f" Max correlation: {np.max(all_correlations):.6f}") print(f" Mean absolute diff (avg): {np.mean(all_mean_abs_diffs):.6f}") + + # Calculate average relative error if threshold is available + if threshold_value and threshold_value > 0: + mean_rel_error = np.mean(all_mean_abs_diffs) / threshold_value + mean_rel_error_ppm = mean_rel_error * 1e6 + print(f" Relative error (avg, scientific): {mean_rel_error:.6e}") + print(f" Relative error (avg, ppm): {mean_rel_error_ppm:.2f} ppm") + print(f" PSNR (avg): {np.mean(all_psnrs):.2f} dB") # Print the directory where files are saved From dde8ba127fdb5ef242afb5180a8a96b10b390bd7 Mon Sep 17 00:00:00 2001 From: "copilot-swe-agent[bot]" <198982749+Copilot@users.noreply.github.com> Date: Thu, 22 May 2025 21:25:36 +0000 Subject: [PATCH 08/21] Fix threshold extraction regex pattern in template matching test Co-authored-by: bHimes <3077528+bHimes@users.noreply.github.com> --- .../match_template/test_template_reproducibility.py | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) diff --git a/scripts/testing/programs/match_template/test_template_reproducibility.py b/scripts/testing/programs/match_template/test_template_reproducibility.py index dad7a8630..a8e752ac2 100755 --- a/scripts/testing/programs/match_template/test_template_reproducibility.py +++ b/scripts/testing/programs/match_template/test_template_reproducibility.py @@ -38,7 +38,7 @@ def extract_threshold_value(hist_file_path): Extract the threshold value from the histogram text file. The threshold value is in the first line of the file, which starts with - "# Expected". This function extracts the numerical value following this marker. + "# Expected threshold = ". This function extracts the numerical value following this marker. Args: hist_file_path (str): Path to the histogram text file @@ -57,8 +57,8 @@ def extract_threshold_value(hist_file_path): first_line = f.readline().strip() # Use regex to extract the threshold value from the line - # Looking for a pattern like "# Expected max is: 0.12345" - match = re.search(r'# Expected max is: ([\d.e+-]+)', first_line) + # Looking for a pattern like "# Expected threshold = 6.90" + match = re.search(r'# Expected threshold = ([\d.e+-]+)', first_line) if not match: raise ValueError(f"Could not find threshold value in {hist_file_path}") From 327a37abb59b1057dae84475002934160e9665a2 Mon Sep 17 00:00:00 2001 From: "copilot-swe-agent[bot]" <198982749+Copilot@users.noreply.github.com> Date: Thu, 22 May 2025 21:37:32 +0000 Subject: [PATCH 09/21] Address PR feedback: improve threshold printing, temp dir management, and CLI options Co-authored-by: bHimes <3077528+bHimes@users.noreply.github.com> --- .../programs/cistem_test_utils/args.py | 76 +++++++++------- .../match_template/.cistem_temp_dirs.log | 1 + .../test_template_reproducibility.py | 90 ++++++++++--------- 3 files changed, 93 insertions(+), 74 deletions(-) create mode 100644 scripts/testing/programs/match_template/.cistem_temp_dirs.log diff --git a/scripts/testing/programs/cistem_test_utils/args.py b/scripts/testing/programs/cistem_test_utils/args.py index cad2fc731..af8049d36 100644 --- a/scripts/testing/programs/cistem_test_utils/args.py +++ b/scripts/testing/programs/cistem_test_utils/args.py @@ -77,7 +77,7 @@ def parse_TM_args(wanted_binary_name): # To help check, we'll keep this dict to use as a check. args_to_check = [] parser = argparse.ArgumentParser(description='Test the k3 rotation binary') - + # Temp directory management arguments temp_management_group = parser.add_argument_group('Temporary Directory Management') temp_management_group.add_argument('--list-temp-dirs', action='store_true', @@ -86,29 +86,29 @@ def parse_TM_args(wanted_binary_name): help='Remove a specific temporary directory by index') temp_management_group.add_argument('--rm-all-temp-dirs', action='store_true', help='Remove all tracked temporary directories') - - # Required argument + + # Binary path argument (required for running tests, optional for temp dir management) parser.add_argument( - '--binary_path', help='Path to the directory with the binary to be tested (Required)', required=True) + '--binary-path', dest='binary_path', help='Path to the directory with the binary to be tested (Required for running tests)', required=False) args_to_check.append('binary_path') - parser.add_argument('--test_data_path', + parser.add_argument('--test-data-path', dest='test_data_path', help='Path to the test data directory (Optional - defaults to /cisTEMdev/cistem_reference_images/TM_tests, then pwd)') args_to_check.append('test_data_path') # Argument for the output file path and prefix default to /tmp - parser.add_argument('--output_file_prefix', + parser.add_argument('--output-file-prefix', dest='output_file_prefix', help='Path and prefix for the output files (Optional - defaults to /tmp)', default='/tmp') args_to_check.append('output_file_prefix') - parser.add_argument('--gpu_idx', default=0, + parser.add_argument('--gpu-idx', dest='gpu_idx', default=0, help='GPU index to use (default: 0)') args_to_check.append('gpu_idx') # add another optional flag to specify that we are using an older version of cisTEM # TODO: for now, just trying to catch the case where we use match_template not match_template_gpu, however, # there could be other cases where we need to be more specific if the input options change more over time. - parser.add_argument('--old_cistem', action='store_true', + parser.add_argument('--old-cistem', dest='old_cistem', action='store_true', help='Use this flag if you are using an older version of cisTEM') args_to_check.append('old_cistem') @@ -117,12 +117,19 @@ def parse_TM_args(wanted_binary_name): help='Use this flag if you are using the cpu version of cisTEM') args_to_check.append('cpu') - parser.add_argument('--fast_fft', action='store_true', default=True, + parser.add_argument('--fast-fft', dest='fast_fft', action='store_true', default=True, help='Use FastFFT implementation (default: True)') args_to_check.append('fast_fft') args = parser.parse_args() + # Check if any temp directory management options are being used + using_temp_management = args.list_temp_dirs or args.rm_temp_dir is not None or args.rm_all_temp_dirs + + # If not using temp management, binary_path is required + if not using_temp_management and not args.binary_path: + parser.error("--binary-path is required when not using temporary directory management options") + args.binary_name = wanted_binary_name args_to_check.append('binary_name') @@ -134,31 +141,40 @@ def parse_TM_args(wanted_binary_name): if not (args.old_cistem or args.cpu): args.binary_name += '_gpu' - # Check if the binary exists - if not os.path.isfile(os.path.join(args.binary_path, args.binary_name)): - print('The binary ' + os.path.join(args.binary_path, - args.binary_name) + ' does not exist') - sys.exit(1) + # Check if any temp directory management options are being used + using_temp_management = args.list_temp_dirs or args.rm_temp_dir is not None or args.rm_all_temp_dirs - # Check if make_template_result binary exists - if not os.path.isfile(os.path.join(args.binary_path, args.results_binary_name)): - print('The binary ' + os.path.join(args.binary_path, - args.results_binary_name) + ' does not exist') - sys.exit(1) + # If not using temp management, binary_path is required and binaries should exist + if not using_temp_management and not args.binary_path: + parser.error("--binary-path is required when not using temporary directory management options") - # if the optional data path is not given, use the default - if args.test_data_path is None: - args.test_data_path = default_data_dir + # Only check for binaries if we're not just managing temp directories + if not using_temp_management: + # Check if the binary exists + if not os.path.isfile(os.path.join(args.binary_path, args.binary_name)): + print('The binary ' + os.path.join(args.binary_path, + args.binary_name) + ' does not exist') + sys.exit(1) - # Check if the test data directory exists - if not os.path.isdir(args.test_data_path): - print('The test data directory [' + - args.test_data_path + '] does not exist') - print('Please provide a valid path to the test data directory as a second argument') - sys.exit(1) + # Check if make_template_result binary exists + if not os.path.isfile(os.path.join(args.binary_path, args.results_binary_name)): + print('The binary ' + os.path.join(args.binary_path, + args.results_binary_name) + ' does not exist') + sys.exit(1) + + # if the optional data path is not given, use the default + if args.test_data_path is None: + args.test_data_path = default_data_dir + + # Check if the test data directory exists + if not os.path.isdir(args.test_data_path): + print('The test data directory [' + + args.test_data_path + '] does not exist') + print('Please provide a valid path to the test data directory as a second argument') + sys.exit(1) - # Check that the wanted output path exists and if not try to make it, if not error - os.makedirs(args.output_file_prefix, exist_ok=True) + # Check that the wanted output path exists and if not try to make it, if not error + os.makedirs(args.output_file_prefix, exist_ok=True) args.args_to_check = args_to_check diff --git a/scripts/testing/programs/match_template/.cistem_temp_dirs.log b/scripts/testing/programs/match_template/.cistem_temp_dirs.log new file mode 100644 index 000000000..0637a088a --- /dev/null +++ b/scripts/testing/programs/match_template/.cistem_temp_dirs.log @@ -0,0 +1 @@ +[] \ No newline at end of file diff --git a/scripts/testing/programs/match_template/test_template_reproducibility.py b/scripts/testing/programs/match_template/test_template_reproducibility.py index a8e752ac2..c071aa397 100755 --- a/scripts/testing/programs/match_template/test_template_reproducibility.py +++ b/scripts/testing/programs/match_template/test_template_reproducibility.py @@ -36,32 +36,32 @@ def extract_threshold_value(hist_file_path): """ Extract the threshold value from the histogram text file. - + The threshold value is in the first line of the file, which starts with "# Expected threshold = ". This function extracts the numerical value following this marker. - + Args: hist_file_path (str): Path to the histogram text file - + Returns: float: The extracted threshold value - + Raises: FileNotFoundError: If the file doesn't exist ValueError: If the threshold value cannot be found or parsed """ if not exists(hist_file_path): raise FileNotFoundError(f"Histogram file not found at {hist_file_path}") - + with open(hist_file_path, 'r') as f: first_line = f.readline().strip() - + # Use regex to extract the threshold value from the line # Looking for a pattern like "# Expected threshold = 6.90" match = re.search(r'# Expected threshold = ([\d.e+-]+)', first_line) if not match: raise ValueError(f"Could not find threshold value in {hist_file_path}") - + threshold_value = float(match.group(1)) return threshold_value @@ -70,17 +70,17 @@ def main(): try: # Parse command-line arguments args = tmArgs.parse_TM_args(wanted_binary_name) - + # Handle temp directory management options if args.list_temp_dirs: temp_dir_manager.print_temp_dirs() return 0 - + if args.rm_temp_dir is not None: success, message = temp_dir_manager.remove_temp_dir(args.rm_temp_dir) print(message) return 0 if success else 1 - + if args.rm_all_temp_dirs: success_count, failure_count = temp_dir_manager.remove_all_temp_dirs() print(f"Successfully removed {success_count} temporary directories.") @@ -91,58 +91,61 @@ def main(): # Create a temporary directory to store our replicate MIPs and track it temp_dir = temp_dir_manager.create_temp_dir(prefix="template_match_reproducibility_") print(f"Temporary directory created at: {temp_dir}") - + # We'll run template matching 3 times to generate replicates elapsed_time = [0, 0, 0] mip_filenames = [] hist_filenames = [] threshold_values = [] - + # Run the template matching 3 times for replicate in range(0, 3): try: # Use Apoferritin dataset with image 0 config = tmArgs.get_config(args, 'Apoferritin', 0, 0) - + # Create a unique output file prefix for each replicate original_prefix = config['output_file_prefix'] config['output_file_prefix'] = join(temp_dir, f"replicate_{replicate+1}") makedirs(config['output_file_prefix'], exist_ok=True) - + # Run the template matching tmp_filename_match_template, tmp_filename_make_template_results = mktmp.make_tmp_runfile(config) elapsed_time[replicate] = runner.run_job(tmp_filename_match_template) runner.run_job(tmp_filename_make_template_results) - + # The MIP file is already in our temp directory with the unique replicate prefix mip_file = join(config['output_file_prefix'], 'mip.mrc') - + # Check if MIP file exists if not exists(mip_file): raise FileNotFoundError(f"MIP file not found at {mip_file}") - + # Get the histogram file path and extract the threshold value hist_file = join(config['output_file_prefix'], 'hist.txt') threshold_value = extract_threshold_value(hist_file) - + mip_filenames.append(mip_file) hist_filenames.append(hist_file) threshold_values.append(threshold_value) - + + # Print threshold value only for the first replicate + if replicate == 0: + print(f"Threshold value: {threshold_value:.6e}") + print(f"Completed replicate {replicate+1}, time: {elapsed_time[replicate]:.2f}s") - print(f" Threshold value: {threshold_value:.6e}") - + except Exception as e: print(f"Error during replicate {replicate+1}: {str(e)}") continue - + # Check if we have all three replicates if len(mip_filenames) != 3: print(f"Warning: Only {len(mip_filenames)} replicates were successfully processed") if len(mip_filenames) < 2: print("At least 2 replicates are required for comparison. Exiting.") return 1 - + # Now load the MIP files and compare them print("\nLoading MIP files for analysis...") mip_data = [] @@ -153,47 +156,46 @@ def main(): print(f"Loaded {filename}, shape: {mrc.data.shape}, dtype: {mrc.data.dtype}") except Exception as e: print(f"Error loading {filename}: {str(e)}") - + # Check if we have at least 2 MIP files loaded if len(mip_data) < 2: print("Failed to load at least 2 MIP files. Cannot perform comparison.") return 1 - + # Verify all threshold values are the same if threshold_values: if not all(abs(v - threshold_values[0]) < 1e-10 for v in threshold_values): print("Warning: Threshold values differ between replicates:") for i, val in enumerate(threshold_values): print(f" Replicate {i+1}: {val:.6e}") - + # Use the first threshold value for calculations threshold_value = threshold_values[0] - print(f"\nThreshold value: {threshold_value:.6e}") else: print("Warning: No threshold values were extracted.") threshold_value = None - + # Compute similarity metrics between replicates print("\nReproducibility Analysis:") print("========================") - + # Pairwise comparisons of available replicates pairs = [(i, j) for i in range(len(mip_data)) for j in range(i+1, len(mip_data))] - + all_correlations = [] all_mean_abs_diffs = [] all_psnrs = [] - + for i, j in pairs: try: # Calculate correlation coefficient correlation = np.corrcoef(mip_data[i].flatten(), mip_data[j].flatten())[0,1] all_correlations.append(correlation) - + # Calculate mean absolute difference mean_abs_diff = np.mean(np.abs(mip_data[i] - mip_data[j])) all_mean_abs_diffs.append(mean_abs_diff) - + # Calculate relative error if threshold value is available if threshold_value and threshold_value > 0: relative_error = mean_abs_diff / threshold_value @@ -201,28 +203,28 @@ def main(): else: relative_error = None relative_error_ppm = None - + # Calculate peak signal to noise ratio (PSNR) mse = np.mean((mip_data[i] - mip_data[j]) ** 2) max_pixel = max(np.max(mip_data[i]), np.max(mip_data[j])) psnr = 20 * np.log10(max_pixel) - 10 * np.log10(mse) if mse > 0 else float('inf') all_psnrs.append(psnr) - + # Print results for this pair print(f"Comparing replicate {i+1} vs {j+1}:") print(f" Correlation coefficient: {correlation:.6f}") print(f" Mean absolute difference: {mean_abs_diff:.6f}") - + # Print relative error if available if relative_error is not None: print(f" Relative error (scientific): {relative_error:.6e}") print(f" Relative error (ppm): {relative_error_ppm:.2f} ppm") - + print(f" Peak signal-to-noise ratio (PSNR): {psnr:.2f} dB") - + except Exception as e: print(f"Error comparing replicates {i+1} and {j+1}: {str(e)}") - + # Overall similarity across all replicates if all_correlations: print("\nOverall reproducibility:") @@ -230,24 +232,24 @@ def main(): print(f" Min correlation: {np.min(all_correlations):.6f}") print(f" Max correlation: {np.max(all_correlations):.6f}") print(f" Mean absolute diff (avg): {np.mean(all_mean_abs_diffs):.6f}") - + # Calculate average relative error if threshold is available if threshold_value and threshold_value > 0: mean_rel_error = np.mean(all_mean_abs_diffs) / threshold_value mean_rel_error_ppm = mean_rel_error * 1e6 print(f" Relative error (avg, scientific): {mean_rel_error:.6e}") print(f" Relative error (avg, ppm): {mean_rel_error_ppm:.2f} ppm") - + print(f" PSNR (avg): {np.mean(all_psnrs):.2f} dB") - + # Print the directory where files are saved print(f"\nMIP files saved in: {temp_dir}") print("To list temp directories: --list-temp-dirs") print("To remove this directory: --rm-temp-dir INDEX") print("To remove all temp directories: --rm-all-temp-dirs") - + return 0 - + except KeyboardInterrupt: print("\nTest interrupted by user.") return 130 From 0556a95fe1d731e085b12c5d5ddf9489d3e88b88 Mon Sep 17 00:00:00 2001 From: "copilot-swe-agent[bot]" <198982749+Copilot@users.noreply.github.com> Date: Fri, 23 May 2025 15:58:59 +0000 Subject: [PATCH 10/21] Implement PR review comments for template matching test Co-authored-by: bHimes <3077528+bHimes@users.noreply.github.com> --- .../programs/cistem_test_utils/args.py | 12 +- .../cistem_test_utils/temp_dir_manager.py | 321 ++++++++++-------- .../cistem_test_utils/threshold_utils.py | 41 +++ .../test_template_reproducibility.py | 100 ++---- 4 files changed, 257 insertions(+), 217 deletions(-) create mode 100644 scripts/testing/programs/cistem_test_utils/threshold_utils.py diff --git a/scripts/testing/programs/cistem_test_utils/args.py b/scripts/testing/programs/cistem_test_utils/args.py index af8049d36..da0b39052 100644 --- a/scripts/testing/programs/cistem_test_utils/args.py +++ b/scripts/testing/programs/cistem_test_utils/args.py @@ -3,6 +3,7 @@ import argparse import toml # If we are in the container +from cistem_test_utils.temp_dir_manager import TempDirManager default_data_dir = '/cisTEMdev/cistem_reference_images/TM_tests' @@ -78,14 +79,9 @@ def parse_TM_args(wanted_binary_name): args_to_check = [] parser = argparse.ArgumentParser(description='Test the k3 rotation binary') - # Temp directory management arguments - temp_management_group = parser.add_argument_group('Temporary Directory Management') - temp_management_group.add_argument('--list-temp-dirs', action='store_true', - help='List all tracked temporary directories') - temp_management_group.add_argument('--rm-temp-dir', type=int, metavar='INDEX', - help='Remove a specific temporary directory by index') - temp_management_group.add_argument('--rm-all-temp-dirs', action='store_true', - help='Remove all tracked temporary directories') + # Add temp directory management arguments using the TempDirManager class + temp_manager = TempDirManager() + temp_manager.add_arguments(parser) # Binary path argument (required for running tests, optional for temp dir management) parser.add_argument( diff --git a/scripts/testing/programs/cistem_test_utils/temp_dir_manager.py b/scripts/testing/programs/cistem_test_utils/temp_dir_manager.py index eb6c69c5d..379ab01cd 100644 --- a/scripts/testing/programs/cistem_test_utils/temp_dir_manager.py +++ b/scripts/testing/programs/cistem_test_utils/temp_dir_manager.py @@ -12,161 +12,206 @@ import datetime import json import shutil +import argparse from typing import List, Dict, Optional, Tuple, Union -# The log file will be a hidden file in the current working directory -LOG_FILE = os.path.join(os.getcwd(), '.cistem_temp_dirs.log') -def create_temp_dir(prefix: str = "cistem_test_", dir: str = "/tmp") -> str: - """ - Create a temporary directory and log its information. - - Args: - prefix: Prefix for the temporary directory name - dir: Parent directory where the temp directory will be created - - Returns: - Path to the created temporary directory - """ - # Create the temporary directory - temp_dir = tempfile.mkdtemp(dir=dir, prefix=prefix) +class TempDirManager: + """Class for managing temporary directories used in cisTEM tests.""" + + def __init__(self): + """Initialize the TempDirManager with a log file in the current working directory.""" + # The log file will be a hidden file in the current working directory + self.log_file = os.path.join(os.getcwd(), '.cistem_temp_dirs.log') + + # Add temporary directory management arguments to an ArgumentParser + self.temp_management_group = None + + def add_arguments(self, parser): + """ + Add temporary directory management arguments to an ArgumentParser. + + Args: + parser: argparse.ArgumentParser instance + """ + self.temp_management_group = parser.add_argument_group('Temporary Directory Management') + self.temp_management_group.add_argument('--list-temp-dirs', action='store_true', + help='List all tracked temporary directories') + self.temp_management_group.add_argument('--rm-temp-dir', type=int, metavar='INDEX', + help='Remove a specific temporary directory by index') + self.temp_management_group.add_argument('--rm-all-temp-dirs', action='store_true', + help='Remove all tracked temporary directories') + + def create_temp_dir(self, prefix: str = "cistem_test_", dir: str = "/tmp") -> str: + """ + Create a temporary directory and log its information. + + Args: + prefix: Prefix for the temporary directory name + dir: Parent directory where the temp directory will be created + + Returns: + Path to the created temporary directory + """ + # Create the temporary directory + temp_dir = tempfile.mkdtemp(dir=dir, prefix=prefix) + + # Log the created directory + log_entry = { + 'path': temp_dir, + 'created_at': datetime.datetime.now().isoformat(), + 'key': os.path.basename(temp_dir) # Use the basename as the unique key + } + + # Read existing log if it exists + log_entries = [] + if os.path.exists(self.log_file): + try: + with open(self.log_file, 'r') as f: + log_entries = json.load(f) + except json.JSONDecodeError: + # If the file is corrupted, start with an empty list + log_entries = [] + + # Add the new entry and write back to the log file + log_entries.append(log_entry) + with open(self.log_file, 'w') as f: + json.dump(log_entries, f, indent=2) + + return temp_dir - # Log the created directory - log_entry = { - 'path': temp_dir, - 'created_at': datetime.datetime.now().isoformat(), - 'key': os.path.basename(temp_dir) # Use the basename as the unique key - } + def list_temp_dirs(self) -> List[Dict]: + """ + List all logged temporary directories. - # Read existing log if it exists - log_entries = [] - if os.path.exists(LOG_FILE): + Returns: + List of dictionaries containing information about each temp directory + """ + if not os.path.exists(self.log_file): + return [] + try: - with open(LOG_FILE, 'r') as f: + with open(self.log_file, 'r') as f: log_entries = json.load(f) + + # Filter out directories that no longer exist + valid_entries = [] + for entry in log_entries: + if os.path.exists(entry['path']): + valid_entries.append(entry) + + # Update the log file if some directories were already removed + if len(valid_entries) != len(log_entries): + with open(self.log_file, 'w') as f: + json.dump(valid_entries, f, indent=2) + + return valid_entries except json.JSONDecodeError: - # If the file is corrupted, start with an empty list - log_entries = [] + # If the file is corrupted, return an empty list + return [] - # Add the new entry and write back to the log file - log_entries.append(log_entry) - with open(LOG_FILE, 'w') as f: - json.dump(log_entries, f, indent=2) + def remove_temp_dir(self, index: int) -> Tuple[bool, str]: + """ + Remove a temporary directory by its index in the log. - return temp_dir - -def list_temp_dirs() -> List[Dict]: - """ - List all logged temporary directories. - - Returns: - List of dictionaries containing information about each temp directory - """ - if not os.path.exists(LOG_FILE): - return [] + Args: + index: The index of the directory to remove (0-based) - try: - with open(LOG_FILE, 'r') as f: - log_entries = json.load(f) + Returns: + Tuple of (success, message) + """ + log_entries = self.list_temp_dirs() + + if not log_entries: + return False, "No temporary directories found in the log." + + if index < 0 or index >= len(log_entries): + return False, f"Invalid index: {index}. Valid range is 0-{len(log_entries)-1}." + + entry = log_entries[index] + path = entry['path'] + + # Check if the directory exists + if not os.path.exists(path): + # Remove from log and return + log_entries.pop(index) + with open(self.log_file, 'w') as f: + json.dump(log_entries, f, indent=2) + return False, f"Directory {path} no longer exists. Removed from log." + + # Remove the directory + try: + shutil.rmtree(path) + # Update the log + log_entries.pop(index) + with open(self.log_file, 'w') as f: + json.dump(log_entries, f, indent=2) + return True, f"Successfully removed {path}" + except Exception as e: + return False, f"Failed to remove {path}: {str(e)}" + + def remove_all_temp_dirs(self) -> Tuple[int, int]: + """ + Remove all temporary directories in the log. + + Returns: + Tuple of (success_count, failure_count) + """ + log_entries = self.list_temp_dirs() + success_count = 0 + failure_count = 0 - # Filter out directories that no longer exist - valid_entries = [] for entry in log_entries: - if os.path.exists(entry['path']): - valid_entries.append(entry) - - # Update the log file if some directories were already removed - if len(valid_entries) != len(log_entries): - with open(LOG_FILE, 'w') as f: - json.dump(valid_entries, f, indent=2) - - return valid_entries - except json.JSONDecodeError: - # If the file is corrupted, return an empty list - return [] + path = entry['path'] + if os.path.exists(path): + try: + shutil.rmtree(path) + success_count += 1 + except: + failure_count += 1 + + # Clear the log file + with open(self.log_file, 'w') as f: + json.dump([], f) + + return success_count, failure_count + + def print_temp_dirs(self) -> None: + """ + Print information about all logged temporary directories. + """ + log_entries = self.list_temp_dirs() + + if not log_entries: + print("No temporary directories found.") + return + + print(f"Found {len(log_entries)} temporary directories:") + print("-" * 80) + for i, entry in enumerate(log_entries): + created_at = datetime.datetime.fromisoformat(entry['created_at']) + formatted_time = created_at.strftime("%Y-%m-%d %H:%M:%S") + print(f"[{i}] {entry['path']}") + print(f" Created: {formatted_time}") + print(f" Key: {entry['key']}") + print("-" * 80) -def remove_temp_dir(index: int) -> Tuple[bool, str]: - """ - Remove a temporary directory by its index in the log. - Args: - index: The index of the directory to remove (0-based) +# Create a singleton instance for backwards compatibility with existing code +_manager = TempDirManager() - Returns: - Tuple of (success, message) - """ - log_entries = list_temp_dirs() - - if not log_entries: - return False, "No temporary directories found in the log." - - if index < 0 or index >= len(log_entries): - return False, f"Invalid index: {index}. Valid range is 0-{len(log_entries)-1}." - - entry = log_entries[index] - path = entry['path'] - - # Check if the directory exists - if not os.path.exists(path): - # Remove from log and return - log_entries.pop(index) - with open(LOG_FILE, 'w') as f: - json.dump(log_entries, f, indent=2) - return False, f"Directory {path} no longer exists. Removed from log." - - # Remove the directory - try: - shutil.rmtree(path) - # Update the log - log_entries.pop(index) - with open(LOG_FILE, 'w') as f: - json.dump(log_entries, f, indent=2) - return True, f"Successfully removed {path}" - except Exception as e: - return False, f"Failed to remove {path}: {str(e)}" +# Provide functions that delegate to the instance methods for backwards compatibility +def create_temp_dir(prefix: str = "cistem_test_", dir: str = "/tmp") -> str: + return _manager.create_temp_dir(prefix, dir) -def remove_all_temp_dirs() -> Tuple[int, int]: - """ - Remove all temporary directories in the log. +def list_temp_dirs() -> List[Dict]: + return _manager.list_temp_dirs() - Returns: - Tuple of (success_count, failure_count) - """ - log_entries = list_temp_dirs() - success_count = 0 - failure_count = 0 - - for entry in log_entries: - path = entry['path'] - if os.path.exists(path): - try: - shutil.rmtree(path) - success_count += 1 - except: - failure_count += 1 - - # Clear the log file - with open(LOG_FILE, 'w') as f: - json.dump([], f) - - return success_count, failure_count +def remove_temp_dir(index: int) -> Tuple[bool, str]: + return _manager.remove_temp_dir(index) + +def remove_all_temp_dirs() -> Tuple[int, int]: + return _manager.remove_all_temp_dirs() def print_temp_dirs() -> None: - """ - Print information about all logged temporary directories. - """ - log_entries = list_temp_dirs() - - if not log_entries: - print("No temporary directories found.") - return - - print(f"Found {len(log_entries)} temporary directories:") - print("-" * 80) - for i, entry in enumerate(log_entries): - created_at = datetime.datetime.fromisoformat(entry['created_at']) - formatted_time = created_at.strftime("%Y-%m-%d %H:%M:%S") - print(f"[{i}] {entry['path']}") - print(f" Created: {formatted_time}") - print(f" Key: {entry['key']}") - print("-" * 80) \ No newline at end of file + return _manager.print_temp_dirs() \ No newline at end of file diff --git a/scripts/testing/programs/cistem_test_utils/threshold_utils.py b/scripts/testing/programs/cistem_test_utils/threshold_utils.py new file mode 100644 index 000000000..a9980757a --- /dev/null +++ b/scripts/testing/programs/cistem_test_utils/threshold_utils.py @@ -0,0 +1,41 @@ +""" +Utilities for handling threshold values in cisTEM outputs. + +This module provides functions for extracting threshold values from cisTEM output files. +""" + +import os +import re + + +def extract_threshold_value(hist_file_path): + """ + Extract the threshold value from the histogram text file. + + The threshold value is in the first line of the file, which starts with + "# Expected threshold = ". This function extracts the numerical value following this marker. + + Args: + hist_file_path (str): Path to the histogram text file + + Returns: + float: The extracted threshold value + + Raises: + FileNotFoundError: If the file doesn't exist + ValueError: If the threshold value cannot be found or parsed + """ + if not os.path.exists(hist_file_path): + raise FileNotFoundError(f"Histogram file not found at {hist_file_path}") + + with open(hist_file_path, 'r') as f: + first_line = f.readline().strip() + + # Use regex to extract the threshold value from the line + # Looking for a pattern like "# Expected threshold = 6.90" + match = re.search(r'# Expected threshold = ([\d.e+-]+)', first_line) + if not match: + raise ValueError(f"Could not find threshold value in {hist_file_path}") + + threshold_value = float(match.group(1)) + return threshold_value \ No newline at end of file diff --git a/scripts/testing/programs/match_template/test_template_reproducibility.py b/scripts/testing/programs/match_template/test_template_reproducibility.py index c071aa397..c4b9e3ad1 100755 --- a/scripts/testing/programs/match_template/test_template_reproducibility.py +++ b/scripts/testing/programs/match_template/test_template_reproducibility.py @@ -2,7 +2,7 @@ """ Template Matching Reproducibility Test -This script runs the cisTEM template matching GPU binary three times using the same input data +This script runs the cisTEM template matching GPU binary multiple times using the same input data and parameters, then analyzes the reproducibility by comparing the resulting MIP images. It uses the Apoferritin dataset for testing, and saves the MIP (Maximum Intensity Projection) @@ -19,7 +19,8 @@ import cistem_test_utils.args as tmArgs import cistem_test_utils.make_tmp_runfile as mktmp import cistem_test_utils.run_job as runner -import cistem_test_utils.temp_dir_manager as temp_dir_manager +from cistem_test_utils.temp_dir_manager import TempDirManager +from cistem_test_utils.threshold_utils import extract_threshold_value import mrcfile import numpy as np import tempfile @@ -32,74 +33,47 @@ # or the --cpu flag is used wanted_binary_name = 'match_template' - -def extract_threshold_value(hist_file_path): - """ - Extract the threshold value from the histogram text file. - - The threshold value is in the first line of the file, which starts with - "# Expected threshold = ". This function extracts the numerical value following this marker. - - Args: - hist_file_path (str): Path to the histogram text file - - Returns: - float: The extracted threshold value - - Raises: - FileNotFoundError: If the file doesn't exist - ValueError: If the threshold value cannot be found or parsed - """ - if not exists(hist_file_path): - raise FileNotFoundError(f"Histogram file not found at {hist_file_path}") - - with open(hist_file_path, 'r') as f: - first_line = f.readline().strip() - - # Use regex to extract the threshold value from the line - # Looking for a pattern like "# Expected threshold = 6.90" - match = re.search(r'# Expected threshold = ([\d.e+-]+)', first_line) - if not match: - raise ValueError(f"Could not find threshold value in {hist_file_path}") - - threshold_value = float(match.group(1)) - return threshold_value +# Define the number of replicates to run +NUM_REPLICATES = 4 def main(): try: + # Create a temp_dir_manager instance + temp_manager = TempDirManager() + # Parse command-line arguments args = tmArgs.parse_TM_args(wanted_binary_name) # Handle temp directory management options if args.list_temp_dirs: - temp_dir_manager.print_temp_dirs() + temp_manager.print_temp_dirs() return 0 if args.rm_temp_dir is not None: - success, message = temp_dir_manager.remove_temp_dir(args.rm_temp_dir) + success, message = temp_manager.remove_temp_dir(args.rm_temp_dir) print(message) return 0 if success else 1 if args.rm_all_temp_dirs: - success_count, failure_count = temp_dir_manager.remove_all_temp_dirs() + success_count, failure_count = temp_manager.remove_all_temp_dirs() print(f"Successfully removed {success_count} temporary directories.") if failure_count > 0: print(f"Failed to remove {failure_count} temporary directories.") return 0 if failure_count == 0 else 1 # Create a temporary directory to store our replicate MIPs and track it - temp_dir = temp_dir_manager.create_temp_dir(prefix="template_match_reproducibility_") + temp_dir = temp_manager.create_temp_dir(prefix="template_match_reproducibility_") print(f"Temporary directory created at: {temp_dir}") - # We'll run template matching 3 times to generate replicates - elapsed_time = [0, 0, 0] + # We'll run template matching for the defined number of replicates + elapsed_time = [0] * NUM_REPLICATES mip_filenames = [] hist_filenames = [] threshold_values = [] - # Run the template matching 3 times - for replicate in range(0, 3): + # Run the template matching for each replicate + for replicate in range(NUM_REPLICATES): try: # Use Apoferritin dataset with image 0 config = tmArgs.get_config(args, 'Apoferritin', 0, 0) @@ -139,12 +113,11 @@ def main(): print(f"Error during replicate {replicate+1}: {str(e)}") continue - # Check if we have all three replicates - if len(mip_filenames) != 3: - print(f"Warning: Only {len(mip_filenames)} replicates were successfully processed") - if len(mip_filenames) < 2: - print("At least 2 replicates are required for comparison. Exiting.") - return 1 + # Check if we have at least 2 replicates + if len(mip_filenames) < 2: + print(f"Error: Only {len(mip_filenames)} replicates were successfully processed") + print("At least 2 replicates are required for comparison. Exiting.") + return 1 # Now load the MIP files and compare them print("\nLoading MIP files for analysis...") @@ -172,19 +145,19 @@ def main(): # Use the first threshold value for calculations threshold_value = threshold_values[0] else: - print("Warning: No threshold values were extracted.") - threshold_value = None + print("Error: No threshold values were extracted.") + return 1 # Compute similarity metrics between replicates print("\nReproducibility Analysis:") print("========================") + print(f"Number of replicates analyzed: {len(mip_data)}") # Pairwise comparisons of available replicates pairs = [(i, j) for i in range(len(mip_data)) for j in range(i+1, len(mip_data))] all_correlations = [] all_mean_abs_diffs = [] - all_psnrs = [] for i, j in pairs: try: @@ -198,29 +171,18 @@ def main(): # Calculate relative error if threshold value is available if threshold_value and threshold_value > 0: - relative_error = mean_abs_diff / threshold_value - relative_error_ppm = relative_error * 1e6 # Parts per million + relative_error_ppm = (mean_abs_diff / threshold_value) * 1e6 # Parts per million else: - relative_error = None relative_error_ppm = None - # Calculate peak signal to noise ratio (PSNR) - mse = np.mean((mip_data[i] - mip_data[j]) ** 2) - max_pixel = max(np.max(mip_data[i]), np.max(mip_data[j])) - psnr = 20 * np.log10(max_pixel) - 10 * np.log10(mse) if mse > 0 else float('inf') - all_psnrs.append(psnr) - # Print results for this pair print(f"Comparing replicate {i+1} vs {j+1}:") print(f" Correlation coefficient: {correlation:.6f}") print(f" Mean absolute difference: {mean_abs_diff:.6f}") # Print relative error if available - if relative_error is not None: - print(f" Relative error (scientific): {relative_error:.6e}") - print(f" Relative error (ppm): {relative_error_ppm:.2f} ppm") - - print(f" Peak signal-to-noise ratio (PSNR): {psnr:.2f} dB") + if relative_error_ppm is not None: + print(f" Relative error: {relative_error_ppm:.2f} ppm (relative to threshold value: {threshold_value:.6e})") except Exception as e: print(f"Error comparing replicates {i+1} and {j+1}: {str(e)}") @@ -235,12 +197,8 @@ def main(): # Calculate average relative error if threshold is available if threshold_value and threshold_value > 0: - mean_rel_error = np.mean(all_mean_abs_diffs) / threshold_value - mean_rel_error_ppm = mean_rel_error * 1e6 - print(f" Relative error (avg, scientific): {mean_rel_error:.6e}") - print(f" Relative error (avg, ppm): {mean_rel_error_ppm:.2f} ppm") - - print(f" PSNR (avg): {np.mean(all_psnrs):.2f} dB") + mean_rel_error_ppm = (np.mean(all_mean_abs_diffs) / threshold_value) * 1e6 + print(f" Relative error (avg): {mean_rel_error_ppm:.2f} ppm (relative to threshold value: {threshold_value:.6e})") # Print the directory where files are saved print(f"\nMIP files saved in: {temp_dir}") From 32505097601fb3f7cba185918571e8c217e40d59 Mon Sep 17 00:00:00 2001 From: "copilot-swe-agent[bot]" <198982749+Copilot@users.noreply.github.com> Date: Fri, 23 May 2025 16:34:35 +0000 Subject: [PATCH 11/21] Implement template matching test improvements Co-authored-by: bHimes <3077528+bHimes@users.noreply.github.com> --- .../programs/cistem_test_utils/args.py | 4 +++ .../cistem_test_utils/make_tmp_runfile.py | 7 ++++-- .../test_template_reproducibility.py | 25 ++++++++----------- 3 files changed, 19 insertions(+), 17 deletions(-) diff --git a/scripts/testing/programs/cistem_test_utils/args.py b/scripts/testing/programs/cistem_test_utils/args.py index da0b39052..f9c441d18 100644 --- a/scripts/testing/programs/cistem_test_utils/args.py +++ b/scripts/testing/programs/cistem_test_utils/args.py @@ -116,6 +116,10 @@ def parse_TM_args(wanted_binary_name): parser.add_argument('--fast-fft', dest='fast_fft', action='store_true', default=True, help='Use FastFFT implementation (default: True)') args_to_check.append('fast_fft') + + parser.add_argument('--binning', dest='binning', type=float, default=1.0, + help='Binning factor for image processing (default: 1.0)') + args_to_check.append('binning') args = parser.parse_args() diff --git a/scripts/testing/programs/cistem_test_utils/make_tmp_runfile.py b/scripts/testing/programs/cistem_test_utils/make_tmp_runfile.py index 0e61e1802..4b9defb55 100644 --- a/scripts/testing/programs/cistem_test_utils/make_tmp_runfile.py +++ b/scripts/testing/programs/cistem_test_utils/make_tmp_runfile.py @@ -21,7 +21,10 @@ def match_template(config): else: use_fast_fft = "no" - high_res_limit = 2*config.get('data')[config.get('img_number')].get('pixel_size') + # Calculate high resolution limit with binning factor + high_res_limit = 2 * config.get('data')[config.get('img_number')].get('pixel_size') * config.get('binning', 1.0) + # Ensure high_res_limit is formatted with at least 4 decimal places + high_res_limit_str = f"{high_res_limit:.4f}" # make a string that will then be used to create a temporary file for feading to stdin # TODO: optional output dir @@ -46,7 +49,7 @@ def match_template(config): str(config.get('data')[config.get('img_number')].get('ctf').get('defocus_2')), str(config.get('data')[config.get('img_number')].get('ctf').get('defocus_angle')), str(config.get('data')[config.get('img_number')].get('ctf').get('extra_phase_shift')), - str(high_res_limit), + high_res_limit_str, str(config.get('out_of_plane_angle')), str(config.get('in_plane_angle')), str(config.get('defocus_range')), diff --git a/scripts/testing/programs/match_template/test_template_reproducibility.py b/scripts/testing/programs/match_template/test_template_reproducibility.py index c4b9e3ad1..d064ba0a7 100755 --- a/scripts/testing/programs/match_template/test_template_reproducibility.py +++ b/scripts/testing/programs/match_template/test_template_reproducibility.py @@ -34,7 +34,7 @@ wanted_binary_name = 'match_template' # Define the number of replicates to run -NUM_REPLICATES = 4 +NUM_REPLICATES = 5 def main(): @@ -77,6 +77,10 @@ def main(): try: # Use Apoferritin dataset with image 0 config = tmArgs.get_config(args, 'Apoferritin', 0, 0) + + # Set the angular sampling step to 3 degrees to speed up testing + config['out_of_plane_angle'] = 3.0 + config['in_plane_angle'] = 3.0 # Create a unique output file prefix for each replicate original_prefix = config['output_file_prefix'] @@ -105,9 +109,9 @@ def main(): # Print threshold value only for the first replicate if replicate == 0: - print(f"Threshold value: {threshold_value:.6e}") + print(f"Threshold value: {threshold_value:.3f}") - print(f"Completed replicate {replicate+1}, time: {elapsed_time[replicate]:.2f}s") + print(f"Completed replicate {replicate+1}/{NUM_REPLICATES}, time: {elapsed_time[replicate]:.2f}s") except Exception as e: print(f"Error during replicate {replicate+1}: {str(e)}") @@ -156,15 +160,10 @@ def main(): # Pairwise comparisons of available replicates pairs = [(i, j) for i in range(len(mip_data)) for j in range(i+1, len(mip_data))] - all_correlations = [] all_mean_abs_diffs = [] for i, j in pairs: try: - # Calculate correlation coefficient - correlation = np.corrcoef(mip_data[i].flatten(), mip_data[j].flatten())[0,1] - all_correlations.append(correlation) - # Calculate mean absolute difference mean_abs_diff = np.mean(np.abs(mip_data[i] - mip_data[j])) all_mean_abs_diffs.append(mean_abs_diff) @@ -177,28 +176,24 @@ def main(): # Print results for this pair print(f"Comparing replicate {i+1} vs {j+1}:") - print(f" Correlation coefficient: {correlation:.6f}") print(f" Mean absolute difference: {mean_abs_diff:.6f}") # Print relative error if available if relative_error_ppm is not None: - print(f" Relative error: {relative_error_ppm:.2f} ppm (relative to threshold value: {threshold_value:.6e})") + print(f" Relative error: {relative_error_ppm:.2f} ppm (relative to threshold value: {threshold_value:.3f})") except Exception as e: print(f"Error comparing replicates {i+1} and {j+1}: {str(e)}") # Overall similarity across all replicates - if all_correlations: + if all_mean_abs_diffs: print("\nOverall reproducibility:") - print(f" Mean correlation: {np.mean(all_correlations):.6f}") - print(f" Min correlation: {np.min(all_correlations):.6f}") - print(f" Max correlation: {np.max(all_correlations):.6f}") print(f" Mean absolute diff (avg): {np.mean(all_mean_abs_diffs):.6f}") # Calculate average relative error if threshold is available if threshold_value and threshold_value > 0: mean_rel_error_ppm = (np.mean(all_mean_abs_diffs) / threshold_value) * 1e6 - print(f" Relative error (avg): {mean_rel_error_ppm:.2f} ppm (relative to threshold value: {threshold_value:.6e})") + print(f" Relative error (avg): {mean_rel_error_ppm:.2f} ppm (relative to threshold value: {threshold_value:.3f})") # Print the directory where files are saved print(f"\nMIP files saved in: {temp_dir}") From dcc9c5e11258ffdee9ff4442ce44958b595cc5e4 Mon Sep 17 00:00:00 2001 From: "copilot-swe-agent[bot]" <198982749+Copilot@users.noreply.github.com> Date: Fri, 23 May 2025 17:15:39 +0000 Subject: [PATCH 12/21] Remove binning CLI flag and add as default config parameter Co-authored-by: bHimes <3077528+bHimes@users.noreply.github.com> --- .cistem_temp_dirs.log | 8 +------- scripts/testing/programs/cistem_test_utils/args.py | 5 +---- .../programs/cistem_test_utils/make_tmp_runfile.py | 2 +- 3 files changed, 3 insertions(+), 12 deletions(-) diff --git a/.cistem_temp_dirs.log b/.cistem_temp_dirs.log index 3c6d38ec9..0637a088a 100644 --- a/.cistem_temp_dirs.log +++ b/.cistem_temp_dirs.log @@ -1,7 +1 @@ -[ - { - "path": "/tmp/cistem_test_ndv2z151", - "created_at": "2025-05-22T21:11:37.778761", - "key": "cistem_test_ndv2z151" - } -] \ No newline at end of file +[] \ No newline at end of file diff --git a/scripts/testing/programs/cistem_test_utils/args.py b/scripts/testing/programs/cistem_test_utils/args.py index f9c441d18..942e621b0 100644 --- a/scripts/testing/programs/cistem_test_utils/args.py +++ b/scripts/testing/programs/cistem_test_utils/args.py @@ -51,6 +51,7 @@ def get_config(args, data_dir: str, ref_number: int, img_number: int): config['padding_factor'] = 1.0 config['mask_radius'] = 0 config['max_threads'] = 2 + config['binning'] = 1.0 # some default search args that may be overwritten in a given test make_template_results config['results_mip_to_use'] = 'mip_scaled.mrc' @@ -116,10 +117,6 @@ def parse_TM_args(wanted_binary_name): parser.add_argument('--fast-fft', dest='fast_fft', action='store_true', default=True, help='Use FastFFT implementation (default: True)') args_to_check.append('fast_fft') - - parser.add_argument('--binning', dest='binning', type=float, default=1.0, - help='Binning factor for image processing (default: 1.0)') - args_to_check.append('binning') args = parser.parse_args() diff --git a/scripts/testing/programs/cistem_test_utils/make_tmp_runfile.py b/scripts/testing/programs/cistem_test_utils/make_tmp_runfile.py index 4b9defb55..541ff9554 100644 --- a/scripts/testing/programs/cistem_test_utils/make_tmp_runfile.py +++ b/scripts/testing/programs/cistem_test_utils/make_tmp_runfile.py @@ -22,7 +22,7 @@ def match_template(config): use_fast_fft = "no" # Calculate high resolution limit with binning factor - high_res_limit = 2 * config.get('data')[config.get('img_number')].get('pixel_size') * config.get('binning', 1.0) + high_res_limit = 2 * config.get('data')[config.get('img_number')].get('pixel_size') * config.get('binning') # Ensure high_res_limit is formatted with at least 4 decimal places high_res_limit_str = f"{high_res_limit:.4f}" From 51450c61b9d29f795ddac5fcfe805aa38ea11137 Mon Sep 17 00:00:00 2001 From: "copilot-swe-agent[bot]" <198982749+Copilot@users.noreply.github.com> Date: Fri, 23 May 2025 17:26:27 +0000 Subject: [PATCH 13/21] Improved Python module imports to replace annoying_hack.py Co-authored-by: bHimes <3077528+bHimes@users.noreply.github.com> --- scripts/testing/README.md | 47 +++++++++++++++++++ .../cistem_test_utils.egg-info/PKG-INFO | 4 ++ .../cistem_test_utils.egg-info/SOURCES.txt | 6 +++ .../dependency_links.txt | 1 + .../cistem_test_utils.egg-info/requires.txt | 3 ++ .../cistem_test_utils.egg-info/top_level.txt | 1 + scripts/testing/programs/cistem_path_setup.py | 22 +++++++++ .../cistem_test_utils.egg-info/PKG-INFO | 4 ++ .../cistem_test_utils.egg-info/SOURCES.txt | 16 +++++++ .../dependency_links.txt | 1 + .../cistem_test_utils.egg-info/requires.txt | 3 ++ .../cistem_test_utils.egg-info/top_level.txt | 1 + .../test_template_reproducibility.py | 6 ++- scripts/testing/programs/setup_env.py | 35 ++++++++++++++ scripts/testing/setup.py | 14 ++++++ 15 files changed, 163 insertions(+), 1 deletion(-) create mode 100644 scripts/testing/README.md create mode 100644 scripts/testing/cistem_test_utils.egg-info/PKG-INFO create mode 100644 scripts/testing/cistem_test_utils.egg-info/SOURCES.txt create mode 100644 scripts/testing/cistem_test_utils.egg-info/dependency_links.txt create mode 100644 scripts/testing/cistem_test_utils.egg-info/requires.txt create mode 100644 scripts/testing/cistem_test_utils.egg-info/top_level.txt create mode 100644 scripts/testing/programs/cistem_path_setup.py create mode 100644 scripts/testing/programs/cistem_test_utils.egg-info/PKG-INFO create mode 100644 scripts/testing/programs/cistem_test_utils.egg-info/SOURCES.txt create mode 100644 scripts/testing/programs/cistem_test_utils.egg-info/dependency_links.txt create mode 100644 scripts/testing/programs/cistem_test_utils.egg-info/requires.txt create mode 100644 scripts/testing/programs/cistem_test_utils.egg-info/top_level.txt create mode 100644 scripts/testing/programs/setup_env.py create mode 100644 scripts/testing/setup.py diff --git a/scripts/testing/README.md b/scripts/testing/README.md new file mode 100644 index 000000000..34ee336c6 --- /dev/null +++ b/scripts/testing/README.md @@ -0,0 +1,47 @@ +# cisTEM Testing Tools + +This directory contains test scripts and utilities for testing cisTEM functionality. + +## Setup + +To properly use these test scripts, you should install the `cistem_test_utils` package in development mode: + +```bash +# Navigate to the testing directory +cd /path/to/cisTEM/scripts/testing + +# Install the package in development mode +pip install -e . +``` + +## Required Dependencies + +The test scripts require the following Python packages: + +- mrcfile +- toml +- numpy + +These can be installed via: + +```bash +pip install mrcfile toml numpy +``` + +## Running Tests + +Each test can be run individually from their respective directories. For example: + +```bash +# Run template matching reproducibility test +python /path/to/cisTEM/scripts/testing/programs/match_template/test_template_reproducibility.py --binary-path /path/to/binaries +``` + +## Temporary Directory Management + +Test scripts that create temporary files use a centralized tracking system to make cleanup easier. +You can list and clean up temporary directories using the following options: + +- `--list-temp-dirs`: List all tracked temporary directories +- `--rm-temp-dir INDEX`: Remove a specific temporary directory by index +- `--rm-all-temp-dirs`: Remove all tracked temporary directories \ No newline at end of file diff --git a/scripts/testing/cistem_test_utils.egg-info/PKG-INFO b/scripts/testing/cistem_test_utils.egg-info/PKG-INFO new file mode 100644 index 000000000..8694c5e1d --- /dev/null +++ b/scripts/testing/cistem_test_utils.egg-info/PKG-INFO @@ -0,0 +1,4 @@ +Metadata-Version: 2.1 +Name: cistem-test-utils +Version: 0.1.0 +Summary: Testing utilities for cisTEM diff --git a/scripts/testing/cistem_test_utils.egg-info/SOURCES.txt b/scripts/testing/cistem_test_utils.egg-info/SOURCES.txt new file mode 100644 index 000000000..87e6c43d5 --- /dev/null +++ b/scripts/testing/cistem_test_utils.egg-info/SOURCES.txt @@ -0,0 +1,6 @@ +setup.py +cistem_test_utils.egg-info/PKG-INFO +cistem_test_utils.egg-info/SOURCES.txt +cistem_test_utils.egg-info/dependency_links.txt +cistem_test_utils.egg-info/requires.txt +cistem_test_utils.egg-info/top_level.txt \ No newline at end of file diff --git a/scripts/testing/cistem_test_utils.egg-info/dependency_links.txt b/scripts/testing/cistem_test_utils.egg-info/dependency_links.txt new file mode 100644 index 000000000..8b1378917 --- /dev/null +++ b/scripts/testing/cistem_test_utils.egg-info/dependency_links.txt @@ -0,0 +1 @@ + diff --git a/scripts/testing/cistem_test_utils.egg-info/requires.txt b/scripts/testing/cistem_test_utils.egg-info/requires.txt new file mode 100644 index 000000000..cbc2e3e74 --- /dev/null +++ b/scripts/testing/cistem_test_utils.egg-info/requires.txt @@ -0,0 +1,3 @@ +mrcfile +numpy +toml diff --git a/scripts/testing/cistem_test_utils.egg-info/top_level.txt b/scripts/testing/cistem_test_utils.egg-info/top_level.txt new file mode 100644 index 000000000..8b1378917 --- /dev/null +++ b/scripts/testing/cistem_test_utils.egg-info/top_level.txt @@ -0,0 +1 @@ + diff --git a/scripts/testing/programs/cistem_path_setup.py b/scripts/testing/programs/cistem_path_setup.py new file mode 100644 index 000000000..148e85d6d --- /dev/null +++ b/scripts/testing/programs/cistem_path_setup.py @@ -0,0 +1,22 @@ +""" +Simple module to replace annoying_hack.py for cisTEM test scripts. + +This module should be imported at the start of any test script that needs +access to cistem_test_utils. It adds the programs directory to the Python path. + +Usage: + import cistem_path_setup # Add this as the first import in test scripts + + # Then import cistem_test_utils modules as normal + import cistem_test_utils.args as tmArgs +""" + +import os +import sys + +# Get the absolute path to the programs directory (parent of this file) +cistem_programs_path = os.path.dirname(os.path.abspath(__file__)) + +# Add to Python path if not already there +if cistem_programs_path not in sys.path: + sys.path.append(cistem_programs_path) \ No newline at end of file diff --git a/scripts/testing/programs/cistem_test_utils.egg-info/PKG-INFO b/scripts/testing/programs/cistem_test_utils.egg-info/PKG-INFO new file mode 100644 index 000000000..8694c5e1d --- /dev/null +++ b/scripts/testing/programs/cistem_test_utils.egg-info/PKG-INFO @@ -0,0 +1,4 @@ +Metadata-Version: 2.1 +Name: cistem-test-utils +Version: 0.1.0 +Summary: Testing utilities for cisTEM diff --git a/scripts/testing/programs/cistem_test_utils.egg-info/SOURCES.txt b/scripts/testing/programs/cistem_test_utils.egg-info/SOURCES.txt new file mode 100644 index 000000000..bed479471 --- /dev/null +++ b/scripts/testing/programs/cistem_test_utils.egg-info/SOURCES.txt @@ -0,0 +1,16 @@ +README.md +setup.py +programs/cistem_test_utils/__init__.py +programs/cistem_test_utils/annoying_hack.py +programs/cistem_test_utils/args.py +programs/cistem_test_utils/compare_coordinates.py +programs/cistem_test_utils/make_tmp_runfile.py +programs/cistem_test_utils/re_run_results_on_mip.py +programs/cistem_test_utils/run_job.py +programs/cistem_test_utils/temp_dir_manager.py +programs/cistem_test_utils/threshold_utils.py +programs/cistem_test_utils.egg-info/PKG-INFO +programs/cistem_test_utils.egg-info/SOURCES.txt +programs/cistem_test_utils.egg-info/dependency_links.txt +programs/cistem_test_utils.egg-info/requires.txt +programs/cistem_test_utils.egg-info/top_level.txt \ No newline at end of file diff --git a/scripts/testing/programs/cistem_test_utils.egg-info/dependency_links.txt b/scripts/testing/programs/cistem_test_utils.egg-info/dependency_links.txt new file mode 100644 index 000000000..8b1378917 --- /dev/null +++ b/scripts/testing/programs/cistem_test_utils.egg-info/dependency_links.txt @@ -0,0 +1 @@ + diff --git a/scripts/testing/programs/cistem_test_utils.egg-info/requires.txt b/scripts/testing/programs/cistem_test_utils.egg-info/requires.txt new file mode 100644 index 000000000..cbc2e3e74 --- /dev/null +++ b/scripts/testing/programs/cistem_test_utils.egg-info/requires.txt @@ -0,0 +1,3 @@ +mrcfile +numpy +toml diff --git a/scripts/testing/programs/cistem_test_utils.egg-info/top_level.txt b/scripts/testing/programs/cistem_test_utils.egg-info/top_level.txt new file mode 100644 index 000000000..7421fb8e4 --- /dev/null +++ b/scripts/testing/programs/cistem_test_utils.egg-info/top_level.txt @@ -0,0 +1 @@ +cistem_test_utils diff --git a/scripts/testing/programs/match_template/test_template_reproducibility.py b/scripts/testing/programs/match_template/test_template_reproducibility.py index d064ba0a7..595dadc4c 100755 --- a/scripts/testing/programs/match_template/test_template_reproducibility.py +++ b/scripts/testing/programs/match_template/test_template_reproducibility.py @@ -12,7 +12,11 @@ measure pixel similarity between replicates with various metrics. """ -import annoying_hack +# Add cisTEM programs directory to Python path +import os, sys +parent_dir = os.path.abspath(os.path.join(os.path.dirname(__file__), "..")) +if parent_dir not in sys.path: + sys.path.insert(0, parent_dir) from os.path import join, exists from os import makedirs diff --git a/scripts/testing/programs/setup_env.py b/scripts/testing/programs/setup_env.py new file mode 100644 index 000000000..3729b0df4 --- /dev/null +++ b/scripts/testing/programs/setup_env.py @@ -0,0 +1,35 @@ +#!/usr/bin/env python3 +""" +Environmental setup utility for cisTEM test scripts. + +This script modifies the Python sys.path to include the necessary directories +so that test scripts can find and import cistem_test_utils without requiring +annoying_hack.py or special environment variables. +""" + +import os +import sys + +def setup_cistem_env(): + """ + Add the necessary directories to sys.path for cisTEM testing. + + This function is designed to be imported and called at the beginning of + test scripts to ensure all the required modules are available. + """ + # Get the absolute path to the programs directory + script_dir = os.path.dirname(os.path.abspath(__file__)) + + # Check if the directories are already in sys.path + if script_dir not in sys.path: + sys.path.append(script_dir) + + return True + +if __name__ == "__main__": + # If run directly, just report success + if setup_cistem_env(): + print("cisTEM environment setup successfully.") + else: + print("Failed to set up cisTEM environment.") + sys.exit(1) \ No newline at end of file diff --git a/scripts/testing/setup.py b/scripts/testing/setup.py new file mode 100644 index 000000000..71a80cd82 --- /dev/null +++ b/scripts/testing/setup.py @@ -0,0 +1,14 @@ +from setuptools import setup, find_packages + +setup( + name="cistem_test_utils", + version="0.1.0", + packages=find_packages(where="programs"), + package_dir={"": "programs"}, + install_requires=[ + "mrcfile", + "toml", + "numpy" + ], + description="Testing utilities for cisTEM", +) \ No newline at end of file From e57c12b5f7b462958c5a4bd66aa66696b65f084e Mon Sep 17 00:00:00 2001 From: "copilot-swe-agent[bot]" <198982749+Copilot@users.noreply.github.com> Date: Fri, 23 May 2025 17:38:29 +0000 Subject: [PATCH 14/21] Modify template matching test to run with multiple binning values Co-authored-by: bHimes <3077528+bHimes@users.noreply.github.com> --- .../test_template_reproducibility.py | 280 +++++++++--------- 1 file changed, 148 insertions(+), 132 deletions(-) diff --git a/scripts/testing/programs/match_template/test_template_reproducibility.py b/scripts/testing/programs/match_template/test_template_reproducibility.py index 595dadc4c..0a9a81475 100755 --- a/scripts/testing/programs/match_template/test_template_reproducibility.py +++ b/scripts/testing/programs/match_template/test_template_reproducibility.py @@ -10,13 +10,12 @@ The script loads the MIP images using the mrcfile package and compares them using numpy to measure pixel similarity between replicates with various metrics. + +The test is run at multiple binning values (1.0, 1.2, 1.6) to verify reproducibility across +different binning configurations. """ -# Add cisTEM programs directory to Python path -import os, sys -parent_dir = os.path.abspath(os.path.join(os.path.dirname(__file__), "..")) -if parent_dir not in sys.path: - sys.path.insert(0, parent_dir) +import annoying_hack from os.path import join, exists from os import makedirs @@ -40,6 +39,9 @@ # Define the number of replicates to run NUM_REPLICATES = 5 +# Define the binning values to test +BINNING_VALUES = [1.0, 1.2, 1.6] + def main(): try: @@ -70,134 +72,148 @@ def main(): temp_dir = temp_manager.create_temp_dir(prefix="template_match_reproducibility_") print(f"Temporary directory created at: {temp_dir}") - # We'll run template matching for the defined number of replicates - elapsed_time = [0] * NUM_REPLICATES - mip_filenames = [] - hist_filenames = [] - threshold_values = [] - - # Run the template matching for each replicate - for replicate in range(NUM_REPLICATES): - try: - # Use Apoferritin dataset with image 0 - config = tmArgs.get_config(args, 'Apoferritin', 0, 0) - - # Set the angular sampling step to 3 degrees to speed up testing - config['out_of_plane_angle'] = 3.0 - config['in_plane_angle'] = 3.0 - - # Create a unique output file prefix for each replicate - original_prefix = config['output_file_prefix'] - config['output_file_prefix'] = join(temp_dir, f"replicate_{replicate+1}") - makedirs(config['output_file_prefix'], exist_ok=True) - - # Run the template matching - tmp_filename_match_template, tmp_filename_make_template_results = mktmp.make_tmp_runfile(config) - elapsed_time[replicate] = runner.run_job(tmp_filename_match_template) - runner.run_job(tmp_filename_make_template_results) - - # The MIP file is already in our temp directory with the unique replicate prefix - mip_file = join(config['output_file_prefix'], 'mip.mrc') - - # Check if MIP file exists - if not exists(mip_file): - raise FileNotFoundError(f"MIP file not found at {mip_file}") - - # Get the histogram file path and extract the threshold value - hist_file = join(config['output_file_prefix'], 'hist.txt') - threshold_value = extract_threshold_value(hist_file) - - mip_filenames.append(mip_file) - hist_filenames.append(hist_file) - threshold_values.append(threshold_value) - - # Print threshold value only for the first replicate - if replicate == 0: - print(f"Threshold value: {threshold_value:.3f}") - - print(f"Completed replicate {replicate+1}/{NUM_REPLICATES}, time: {elapsed_time[replicate]:.2f}s") - - except Exception as e: - print(f"Error during replicate {replicate+1}: {str(e)}") - continue - - # Check if we have at least 2 replicates - if len(mip_filenames) < 2: - print(f"Error: Only {len(mip_filenames)} replicates were successfully processed") - print("At least 2 replicates are required for comparison. Exiting.") - return 1 - - # Now load the MIP files and compare them - print("\nLoading MIP files for analysis...") - mip_data = [] - for filename in mip_filenames: - try: - with mrcfile.open(filename) as mrc: - mip_data.append(mrc.data) - print(f"Loaded {filename}, shape: {mrc.data.shape}, dtype: {mrc.data.dtype}") - except Exception as e: - print(f"Error loading {filename}: {str(e)}") - - # Check if we have at least 2 MIP files loaded - if len(mip_data) < 2: - print("Failed to load at least 2 MIP files. Cannot perform comparison.") - return 1 - - # Verify all threshold values are the same - if threshold_values: - if not all(abs(v - threshold_values[0]) < 1e-10 for v in threshold_values): - print("Warning: Threshold values differ between replicates:") - for i, val in enumerate(threshold_values): - print(f" Replicate {i+1}: {val:.6e}") - - # Use the first threshold value for calculations - threshold_value = threshold_values[0] - else: - print("Error: No threshold values were extracted.") - return 1 - - # Compute similarity metrics between replicates - print("\nReproducibility Analysis:") - print("========================") - print(f"Number of replicates analyzed: {len(mip_data)}") - - # Pairwise comparisons of available replicates - pairs = [(i, j) for i in range(len(mip_data)) for j in range(i+1, len(mip_data))] - - all_mean_abs_diffs = [] - - for i, j in pairs: - try: - # Calculate mean absolute difference - mean_abs_diff = np.mean(np.abs(mip_data[i] - mip_data[j])) - all_mean_abs_diffs.append(mean_abs_diff) - - # Calculate relative error if threshold value is available + # Test each binning value + for binning_value in BINNING_VALUES: + print(f"\n{'-'*80}") + print(f"Testing with binning value: {binning_value}") + print(f"{'-'*80}") + + # Create a subdirectory for this binning value + binning_dir = join(temp_dir, f"binning_{binning_value}") + os.makedirs(binning_dir, exist_ok=True) + + # We'll run template matching for the defined number of replicates + elapsed_time = [0] * NUM_REPLICATES + mip_filenames = [] + hist_filenames = [] + threshold_values = [] + + # Run the template matching for each replicate + for replicate in range(NUM_REPLICATES): + try: + # Use Apoferritin dataset with image 0 + config = tmArgs.get_config(args, 'Apoferritin', 0, 0) + + # Set the angular sampling step to 3 degrees to speed up testing + config['out_of_plane_angle'] = 3.0 + config['in_plane_angle'] = 3.0 + + # Set the binning value + config['binning'] = binning_value + + # Create a unique output file prefix for each replicate + original_prefix = config['output_file_prefix'] + config['output_file_prefix'] = join(binning_dir, f"replicate_{replicate+1}") + makedirs(config['output_file_prefix'], exist_ok=True) + + # Run the template matching + tmp_filename_match_template, tmp_filename_make_template_results = mktmp.make_tmp_runfile(config) + elapsed_time[replicate] = runner.run_job(tmp_filename_match_template) + runner.run_job(tmp_filename_make_template_results) + + # The MIP file is already in our temp directory with the unique replicate prefix + mip_file = join(config['output_file_prefix'], 'mip.mrc') + + # Check if MIP file exists + if not exists(mip_file): + raise FileNotFoundError(f"MIP file not found at {mip_file}") + + # Get the histogram file path and extract the threshold value + hist_file = join(config['output_file_prefix'], 'hist.txt') + threshold_value = extract_threshold_value(hist_file) + + mip_filenames.append(mip_file) + hist_filenames.append(hist_file) + threshold_values.append(threshold_value) + + # Print threshold value only for the first replicate + if replicate == 0: + print(f"Threshold value: {threshold_value:.3f}") + + print(f"Completed replicate {replicate+1}/{NUM_REPLICATES}, time: {elapsed_time[replicate]:.2f}s") + + except Exception as e: + print(f"Error during replicate {replicate+1}: {str(e)}") + continue + + # Check if we have at least 2 replicates + if len(mip_filenames) < 2: + print(f"Error: Only {len(mip_filenames)} replicates were successfully processed") + print("At least 2 replicates are required for comparison. Exiting.") + return 1 + + # Now load the MIP files and compare them + print("\nLoading MIP files for analysis...") + mip_data = [] + for filename in mip_filenames: + try: + with mrcfile.open(filename) as mrc: + mip_data.append(mrc.data) + print(f"Loaded {filename}, shape: {mrc.data.shape}, dtype: {mrc.data.dtype}") + except Exception as e: + print(f"Error loading {filename}: {str(e)}") + + # Check if we have at least 2 MIP files loaded + if len(mip_data) < 2: + print("Failed to load at least 2 MIP files. Cannot perform comparison.") + return 1 + + # Verify all threshold values are the same + if threshold_values: + if not all(abs(v - threshold_values[0]) < 1e-10 for v in threshold_values): + print("Warning: Threshold values differ between replicates:") + for i, val in enumerate(threshold_values): + print(f" Replicate {i+1}: {val:.6e}") + + # Use the first threshold value for calculations + threshold_value = threshold_values[0] + else: + print("Error: No threshold values were extracted.") + return 1 + + # Compute similarity metrics between replicates + print("\nReproducibility Analysis:") + print("========================") + print(f"Number of replicates analyzed: {len(mip_data)}") + print(f"Binning value: {binning_value}") + + # Pairwise comparisons of available replicates + pairs = [(i, j) for i in range(len(mip_data)) for j in range(i+1, len(mip_data))] + + all_mean_abs_diffs = [] + + for i, j in pairs: + try: + # Calculate mean absolute difference + mean_abs_diff = np.mean(np.abs(mip_data[i] - mip_data[j])) + all_mean_abs_diffs.append(mean_abs_diff) + + # Calculate relative error if threshold value is available + if threshold_value and threshold_value > 0: + relative_error_ppm = (mean_abs_diff / threshold_value) * 1e6 # Parts per million + else: + relative_error_ppm = None + + # Print results for this pair + print(f"Comparing replicate {i+1} vs {j+1}:") + print(f" Mean absolute difference: {mean_abs_diff:.6f}") + + # Print relative error if available + if relative_error_ppm is not None: + print(f" Relative error: {relative_error_ppm:.2f} ppm (relative to threshold value: {threshold_value:.3f})") + + except Exception as e: + print(f"Error comparing replicates {i+1} and {j+1}: {str(e)}") + + # Overall similarity across all replicates + if all_mean_abs_diffs: + print("\nOverall reproducibility:") + print(f" Mean absolute diff (avg): {np.mean(all_mean_abs_diffs):.6f}") + + # Calculate average relative error if threshold is available if threshold_value and threshold_value > 0: - relative_error_ppm = (mean_abs_diff / threshold_value) * 1e6 # Parts per million - else: - relative_error_ppm = None - - # Print results for this pair - print(f"Comparing replicate {i+1} vs {j+1}:") - print(f" Mean absolute difference: {mean_abs_diff:.6f}") - - # Print relative error if available - if relative_error_ppm is not None: - print(f" Relative error: {relative_error_ppm:.2f} ppm (relative to threshold value: {threshold_value:.3f})") - - except Exception as e: - print(f"Error comparing replicates {i+1} and {j+1}: {str(e)}") - - # Overall similarity across all replicates - if all_mean_abs_diffs: - print("\nOverall reproducibility:") - print(f" Mean absolute diff (avg): {np.mean(all_mean_abs_diffs):.6f}") - - # Calculate average relative error if threshold is available - if threshold_value and threshold_value > 0: - mean_rel_error_ppm = (np.mean(all_mean_abs_diffs) / threshold_value) * 1e6 - print(f" Relative error (avg): {mean_rel_error_ppm:.2f} ppm (relative to threshold value: {threshold_value:.3f})") + mean_rel_error_ppm = (np.mean(all_mean_abs_diffs) / threshold_value) * 1e6 + print(f" Relative error (avg): {mean_rel_error_ppm:.2f} ppm (relative to threshold value: {threshold_value:.3f})") # Print the directory where files are saved print(f"\nMIP files saved in: {temp_dir}") From eb0b6a0d8865fa5f28044f3f4cef7ee3452925c0 Mon Sep 17 00:00:00 2001 From: "copilot-swe-agent[bot]" <198982749+Copilot@users.noreply.github.com> Date: Sat, 24 May 2025 11:59:06 +0000 Subject: [PATCH 15/21] Create reusable image replicate analysis tools Co-authored-by: bHimes <3077528+bHimes@users.noreply.github.com> --- .../programs/cistem_test_utils/__init__.py | 2 + .../image_replicate_analysis.py | 182 ++++++++++++++++++ .../programs/image_replicate_analysis.py | 118 ++++++++++++ .../test_template_reproducibility.py | 80 ++------ 4 files changed, 323 insertions(+), 59 deletions(-) create mode 100644 scripts/testing/programs/cistem_test_utils/image_replicate_analysis.py create mode 100755 scripts/testing/programs/image_replicate_analysis.py diff --git a/scripts/testing/programs/cistem_test_utils/__init__.py b/scripts/testing/programs/cistem_test_utils/__init__.py index 756dbcee7..9d35520ec 100644 --- a/scripts/testing/programs/cistem_test_utils/__init__.py +++ b/scripts/testing/programs/cistem_test_utils/__init__.py @@ -4,3 +4,5 @@ from . import make_tmp_runfile from . import run_job from . import temp_dir_manager +from . import threshold_utils +from . import image_replicate_analysis diff --git a/scripts/testing/programs/cistem_test_utils/image_replicate_analysis.py b/scripts/testing/programs/cistem_test_utils/image_replicate_analysis.py new file mode 100644 index 000000000..433ead096 --- /dev/null +++ b/scripts/testing/programs/cistem_test_utils/image_replicate_analysis.py @@ -0,0 +1,182 @@ +""" +Image Replicate Analysis utilities for cisTEM tests. + +This module provides a class for comparing multiple replicate MRC images and calculating +similarity metrics between them. +""" + +import numpy as np +import mrcfile +import os +from typing import List, Dict, Tuple, Optional, Union + + +class ImageReplicateAnalysis: + """Class for analyzing replicate MRC images and calculating similarity metrics.""" + + def __init__(self, image_filenames: List[str], threshold_value: float = None): + """ + Initialize the ImageReplicateAnalysis with a list of image filenames and threshold value. + + Args: + image_filenames: List of MRC image files to analyze + threshold_value: Threshold value for relative error calculations (must be between 0 and 100) + + Raises: + ValueError: If threshold_value is not between 0 and 100 + ValueError: If fewer than 2 image filenames are provided + """ + if len(image_filenames) < 2: + raise ValueError("At least 2 image filenames are required for comparison") + + self.image_filenames = image_filenames + + # Validate threshold value if provided + if threshold_value is not None: + if not isinstance(threshold_value, (int, float)) or threshold_value <= 0 or threshold_value > 100: + raise ValueError("Threshold value must be a positive number between 0 and 100") + + self.threshold_value = threshold_value + self.image_data = [] + self.image_shapes = [] + self.image_dtypes = [] + + def load_images(self) -> bool: + """ + Load all image files and verify they have the same dimensions. + + Returns: + bool: True if all images were loaded successfully with matching dimensions + + Raises: + FileNotFoundError: If any image file cannot be found + ValueError: If image dimensions do not match + """ + self.image_data = [] + self.image_shapes = [] + self.image_dtypes = [] + + # Load all images + for filename in self.image_filenames: + if not os.path.exists(filename): + raise FileNotFoundError(f"Image file not found: {filename}") + + try: + with mrcfile.open(filename) as mrc: + self.image_data.append(mrc.data) + self.image_shapes.append(mrc.data.shape) + self.image_dtypes.append(mrc.data.dtype) + except Exception as e: + raise IOError(f"Error loading {filename}: {str(e)}") + + # Check that all images have the same dimensions + if len(set(str(shape) for shape in self.image_shapes)) > 1: + raise ValueError(f"Image dimensions do not match: {self.image_shapes}") + + return True + + def analyze_replicates(self) -> Dict: + """ + Analyze all replicate images and calculate similarity metrics. + + Returns: + Dict: Dictionary containing pairwise and overall similarity metrics + """ + if not self.image_data: + self.load_images() + + num_replicates = len(self.image_data) + + # Generate all pairwise comparisons + pairs = [(i, j) for i in range(num_replicates) for j in range(i+1, num_replicates)] + + results = { + "num_replicates": num_replicates, + "threshold_value": self.threshold_value, + "pairwise_comparisons": [], + "overall": {} + } + + all_mean_abs_diffs = [] + + # Calculate metrics for each pair of images + for i, j in pairs: + try: + # Calculate mean absolute difference + mean_abs_diff = np.mean(np.abs(self.image_data[i] - self.image_data[j])) + all_mean_abs_diffs.append(mean_abs_diff) + + # Calculate relative error if threshold value is available + if self.threshold_value and self.threshold_value > 0: + relative_error_ppm = (mean_abs_diff / self.threshold_value) * 1e6 # Parts per million + else: + relative_error_ppm = None + + comparison_result = { + "replicate_1": i + 1, + "replicate_2": j + 1, + "mean_abs_diff": mean_abs_diff, + "relative_error_ppm": relative_error_ppm + } + + results["pairwise_comparisons"].append(comparison_result) + + except Exception as e: + print(f"Error comparing replicates {i+1} and {j+1}: {str(e)}") + + # Calculate overall metrics across all comparisons + if all_mean_abs_diffs: + results["overall"]["mean_abs_diff_avg"] = np.mean(all_mean_abs_diffs) + results["overall"]["mean_abs_diff_min"] = np.min(all_mean_abs_diffs) + results["overall"]["mean_abs_diff_max"] = np.max(all_mean_abs_diffs) + + # Calculate average relative error if threshold is available + if self.threshold_value and self.threshold_value > 0: + results["overall"]["relative_error_ppm_avg"] = (np.mean(all_mean_abs_diffs) / self.threshold_value) * 1e6 + results["overall"]["relative_error_ppm_min"] = (np.min(all_mean_abs_diffs) / self.threshold_value) * 1e6 + results["overall"]["relative_error_ppm_max"] = (np.max(all_mean_abs_diffs) / self.threshold_value) * 1e6 + + return results + + def print_analysis(self, results: Optional[Dict] = None) -> None: + """ + Print the replicate analysis results in a formatted way. + + Args: + results: Optional results dictionary from analyze_replicates(). + If None, will run analyze_replicates() internally. + """ + if results is None: + results = self.analyze_replicates() + + num_replicates = results["num_replicates"] + threshold_value = results["threshold_value"] + + print("\nReproducibility Analysis:") + print("========================") + print(f"Number of replicates analyzed: {num_replicates}") + + if threshold_value is not None: + print(f"Threshold value: {threshold_value:.3f}") + + # Print pairwise comparisons + for comparison in results["pairwise_comparisons"]: + i = comparison["replicate_1"] + j = comparison["replicate_2"] + mean_abs_diff = comparison["mean_abs_diff"] + relative_error_ppm = comparison["relative_error_ppm"] + + print(f"\nComparing replicate {i} vs {j}:") + print(f" Mean absolute difference: {mean_abs_diff:.6f}") + + if relative_error_ppm is not None: + print(f" Relative error: {relative_error_ppm:.2f} ppm (relative to threshold value: {threshold_value:.3f})") + + # Print overall metrics + if "overall" in results and results["overall"]: + print("\nOverall reproducibility:") + print(f" Mean absolute diff (avg): {results['overall']['mean_abs_diff_avg']:.6f}") + + if threshold_value is not None: + print(f" Relative error (avg): {results['overall']['relative_error_ppm_avg']:.2f} ppm " + f"(relative to threshold value: {threshold_value:.3f})") \ No newline at end of file diff --git a/scripts/testing/programs/image_replicate_analysis.py b/scripts/testing/programs/image_replicate_analysis.py new file mode 100755 index 000000000..ba5a6728c --- /dev/null +++ b/scripts/testing/programs/image_replicate_analysis.py @@ -0,0 +1,118 @@ +#!/usr/bin/env python3 +""" +Image Replicate Analysis Tool + +This script analyzes multiple replicate MRC images and calculates similarity metrics +between them. It accepts a file containing a list of image filenames and a threshold +value for calculating relative error metrics. + +Example usage: + python image_replicate_analysis.py --image-list /path/to/image_list.txt --threshold 6.90 +""" + +import sys +import os +import argparse + +# Handle Python module import paths +script_dir = os.path.dirname(os.path.abspath(__file__)) +sys.path.append(os.path.dirname(script_dir)) + +try: + from cistem_test_utils.image_replicate_analysis import ImageReplicateAnalysis +except ImportError: + # Fall back to annoying_hack.py approach for compatibility + sys.path.insert(0, os.path.join(script_dir)) + import annoying_hack + from cistem_test_utils.image_replicate_analysis import ImageReplicateAnalysis + + +def parse_args(): + """Parse command line arguments.""" + parser = argparse.ArgumentParser( + description="Analyze replicate MRC images and calculate similarity metrics" + ) + + parser.add_argument( + "--image-list", + required=True, + help="Path to a file containing a list of MRC image filenames (one per line)" + ) + parser.add_argument( + "--threshold", + type=float, + required=True, + help="Threshold value for relative error calculations (must be between 0 and 100)" + ) + + return parser.parse_args() + + +def read_image_list(file_path): + """ + Read image filenames from a text file. + + Args: + file_path: Path to the text file containing image filenames + + Returns: + List of image filenames + + Raises: + FileNotFoundError: If the file doesn't exist + ValueError: If the file is empty or contains invalid entries + """ + if not os.path.exists(file_path): + raise FileNotFoundError(f"Image list file not found: {file_path}") + + with open(file_path, 'r') as f: + image_files = [line.strip() for line in f if line.strip()] + + if not image_files: + raise ValueError(f"No image filenames found in {file_path}") + + # Check that all image files exist + missing_files = [f for f in image_files if not os.path.exists(f)] + if missing_files: + raise FileNotFoundError(f"The following image files were not found: {', '.join(missing_files)}") + + return image_files + + +def main(): + """Main entry point for the script.""" + try: + args = parse_args() + + # Read image filenames from the list file + image_files = read_image_list(args.image_list) + print(f"Found {len(image_files)} image files in {args.image_list}") + + # Create image replicate analyzer instance + analyzer = ImageReplicateAnalysis(image_files, args.threshold) + + # Load images and verify dimensions + print("Loading images and verifying dimensions...") + analyzer.load_images() + + # Run the analysis + print("Analyzing replicate images...") + results = analyzer.analyze_replicates() + + # Print the results + analyzer.print_analysis(results) + + return 0 + + except KeyboardInterrupt: + print("\nAnalysis interrupted by user.") + return 130 + except Exception as e: + print(f"Error: {str(e)}") + import traceback + traceback.print_exc() + return 1 + + +if __name__ == "__main__": + sys.exit(main()) \ No newline at end of file diff --git a/scripts/testing/programs/match_template/test_template_reproducibility.py b/scripts/testing/programs/match_template/test_template_reproducibility.py index 0a9a81475..ef04276b0 100755 --- a/scripts/testing/programs/match_template/test_template_reproducibility.py +++ b/scripts/testing/programs/match_template/test_template_reproducibility.py @@ -24,6 +24,7 @@ import cistem_test_utils.run_job as runner from cistem_test_utils.temp_dir_manager import TempDirManager from cistem_test_utils.threshold_utils import extract_threshold_value +from cistem_test_utils.image_replicate_analysis import ImageReplicateAnalysis import mrcfile import numpy as np import tempfile @@ -142,22 +143,6 @@ def main(): print("At least 2 replicates are required for comparison. Exiting.") return 1 - # Now load the MIP files and compare them - print("\nLoading MIP files for analysis...") - mip_data = [] - for filename in mip_filenames: - try: - with mrcfile.open(filename) as mrc: - mip_data.append(mrc.data) - print(f"Loaded {filename}, shape: {mrc.data.shape}, dtype: {mrc.data.dtype}") - except Exception as e: - print(f"Error loading {filename}: {str(e)}") - - # Check if we have at least 2 MIP files loaded - if len(mip_data) < 2: - print("Failed to load at least 2 MIP files. Cannot perform comparison.") - return 1 - # Verify all threshold values are the same if threshold_values: if not all(abs(v - threshold_values[0]) < 1e-10 for v in threshold_values): @@ -170,50 +155,27 @@ def main(): else: print("Error: No threshold values were extracted.") return 1 - - # Compute similarity metrics between replicates - print("\nReproducibility Analysis:") - print("========================") - print(f"Number of replicates analyzed: {len(mip_data)}") + print(f"Binning value: {binning_value}") - - # Pairwise comparisons of available replicates - pairs = [(i, j) for i in range(len(mip_data)) for j in range(i+1, len(mip_data))] - - all_mean_abs_diffs = [] - - for i, j in pairs: - try: - # Calculate mean absolute difference - mean_abs_diff = np.mean(np.abs(mip_data[i] - mip_data[j])) - all_mean_abs_diffs.append(mean_abs_diff) - - # Calculate relative error if threshold value is available - if threshold_value and threshold_value > 0: - relative_error_ppm = (mean_abs_diff / threshold_value) * 1e6 # Parts per million - else: - relative_error_ppm = None - - # Print results for this pair - print(f"Comparing replicate {i+1} vs {j+1}:") - print(f" Mean absolute difference: {mean_abs_diff:.6f}") - - # Print relative error if available - if relative_error_ppm is not None: - print(f" Relative error: {relative_error_ppm:.2f} ppm (relative to threshold value: {threshold_value:.3f})") - - except Exception as e: - print(f"Error comparing replicates {i+1} and {j+1}: {str(e)}") - - # Overall similarity across all replicates - if all_mean_abs_diffs: - print("\nOverall reproducibility:") - print(f" Mean absolute diff (avg): {np.mean(all_mean_abs_diffs):.6f}") - - # Calculate average relative error if threshold is available - if threshold_value and threshold_value > 0: - mean_rel_error_ppm = (np.mean(all_mean_abs_diffs) / threshold_value) * 1e6 - print(f" Relative error (avg): {mean_rel_error_ppm:.2f} ppm (relative to threshold value: {threshold_value:.3f})") + + # Use the ImageReplicateAnalysis class for image comparison + try: + # Create analyzer instance with MIP filenames and threshold value + analyzer = ImageReplicateAnalysis(mip_filenames, threshold_value) + + # Load the images + print("\nLoading MIP files for analysis...") + analyzer.load_images() + + # Run the analysis and print the results + results = analyzer.analyze_replicates() + analyzer.print_analysis(results) + + except Exception as e: + print(f"Error during image replicate analysis: {str(e)}") + import traceback + traceback.print_exc() + return 1 # Print the directory where files are saved print(f"\nMIP files saved in: {temp_dir}") From 9004eb5656c14adb539867e6072d7a7b4697005a Mon Sep 17 00:00:00 2001 From: "copilot-swe-agent[bot]" <198982749+Copilot@users.noreply.github.com> Date: Sat, 24 May 2025 12:14:46 +0000 Subject: [PATCH 16/21] Modify test_template_reproducibility.py to add fast development mode Co-authored-by: bHimes <3077528+bHimes@users.noreply.github.com> --- .../test_template_reproducibility.py | 62 +++++++++++-------- 1 file changed, 37 insertions(+), 25 deletions(-) diff --git a/scripts/testing/programs/match_template/test_template_reproducibility.py b/scripts/testing/programs/match_template/test_template_reproducibility.py index ef04276b0..edcf3ef25 100755 --- a/scripts/testing/programs/match_template/test_template_reproducibility.py +++ b/scripts/testing/programs/match_template/test_template_reproducibility.py @@ -11,8 +11,12 @@ The script loads the MIP images using the mrcfile package and compares them using numpy to measure pixel similarity between replicates with various metrics. -The test is run at multiple binning values (1.0, 1.2, 1.6) to verify reproducibility across -different binning configurations. +The script has two operational modes controlled by the FAST_DEVELOPMENT_CONDITIONS flag: +1. When True: Uses only binning=1.2 with 5.0 degree angular sampling (faster development mode) +2. When False: Tests all binning values (1.0, 1.2, 1.6) with 3.0 degree angular sampling + +Instead of performing analysis directly, this script outputs commands for the user to run +the image_replicate_analysis.py tool separately. """ import annoying_hack @@ -40,8 +44,14 @@ # Define the number of replicates to run NUM_REPLICATES = 5 -# Define the binning values to test -BINNING_VALUES = [1.0, 1.2, 1.6] +# Define whether to use fast development conditions +FAST_DEVELOPMENT_CONDITIONS = True + +# Define the binning values to test - modified based on development conditions +if FAST_DEVELOPMENT_CONDITIONS: + BINNING_VALUES = [1.2] # Only use 1.2 binning for faster development +else: + BINNING_VALUES = [1.0, 1.2, 1.6] def main(): @@ -95,9 +105,13 @@ def main(): # Use Apoferritin dataset with image 0 config = tmArgs.get_config(args, 'Apoferritin', 0, 0) - # Set the angular sampling step to 3 degrees to speed up testing - config['out_of_plane_angle'] = 3.0 - config['in_plane_angle'] = 3.0 + # Set the angular sampling step based on development conditions + if FAST_DEVELOPMENT_CONDITIONS: + config['out_of_plane_angle'] = 5.0 + config['in_plane_angle'] = 5.0 + else: + config['out_of_plane_angle'] = 3.0 + config['in_plane_angle'] = 3.0 # Set the binning value config['binning'] = binning_value @@ -158,24 +172,22 @@ def main(): print(f"Binning value: {binning_value}") - # Use the ImageReplicateAnalysis class for image comparison - try: - # Create analyzer instance with MIP filenames and threshold value - analyzer = ImageReplicateAnalysis(mip_filenames, threshold_value) - - # Load the images - print("\nLoading MIP files for analysis...") - analyzer.load_images() - - # Run the analysis and print the results - results = analyzer.analyze_replicates() - analyzer.print_analysis(results) - - except Exception as e: - print(f"Error during image replicate analysis: {str(e)}") - import traceback - traceback.print_exc() - return 1 + # Write MIP filenames to a file in the binning directory + mip_list_file = join(binning_dir, "mip_files.txt") + with open(mip_list_file, 'w') as f: + for mip_file in mip_filenames: + f.write(f"{mip_file}\n") + + # Instead of running analysis directly, print commands for the user to run + print(f"\n{'-'*80}") + print(f"Replicate generation complete for binning {binning_value}") + print(f"MIP filenames saved to: {mip_list_file}") + print(f"{'-'*80}") + print("\nTo analyze the replicate results, run:") + print(f"python3 scripts/testing/programs/cistem_test_utils/image_replicate_analysis.py --image-list {mip_list_file} --threshold {threshold_value:.3f}") + print("\nOr to extract the threshold value directly from histogram files:") + print(f"threshold=$(awk '/Threshold/{{print $NF; exit}}' {hist_filenames[0]})") + print(f"python3 scripts/testing/programs/cistem_test_utils/image_replicate_analysis.py --image-list {mip_list_file} --threshold $threshold") # Print the directory where files are saved print(f"\nMIP files saved in: {temp_dir}") From baf1bce568cf958639938eb04a05c51fda520284 Mon Sep 17 00:00:00 2001 From: "copilot-swe-agent[bot]" <198982749+Copilot@users.noreply.github.com> Date: Sat, 24 May 2025 12:27:41 +0000 Subject: [PATCH 17/21] Update template matching test to provide bash commands for analysis Co-authored-by: bHimes <3077528+bHimes@users.noreply.github.com> --- .../test_template_reproducibility.py | 24 +++++++++++-------- 1 file changed, 14 insertions(+), 10 deletions(-) diff --git a/scripts/testing/programs/match_template/test_template_reproducibility.py b/scripts/testing/programs/match_template/test_template_reproducibility.py index edcf3ef25..78e25582a 100755 --- a/scripts/testing/programs/match_template/test_template_reproducibility.py +++ b/scripts/testing/programs/match_template/test_template_reproducibility.py @@ -172,22 +172,26 @@ def main(): print(f"Binning value: {binning_value}") - # Write MIP filenames to a file in the binning directory - mip_list_file = join(binning_dir, "mip_files.txt") - with open(mip_list_file, 'w') as f: - for mip_file in mip_filenames: - f.write(f"{mip_file}\n") - # Instead of running analysis directly, print commands for the user to run print(f"\n{'-'*80}") print(f"Replicate generation complete for binning {binning_value}") - print(f"MIP filenames saved to: {mip_list_file}") print(f"{'-'*80}") - print("\nTo analyze the replicate results, run:") - print(f"python3 scripts/testing/programs/cistem_test_utils/image_replicate_analysis.py --image-list {mip_list_file} --threshold {threshold_value:.3f}") + + # Suggest a list file location but don't create it - let user do it with wildcards + mip_list_file = join(binning_dir, "mip_files.txt") + + print("\nTo create a list of MIP files for analysis:") + print(f"find {binning_dir} -name 'mip.mrc' > {mip_list_file}") + + print("\nOr to analyze other image types, like theta maps:") + print(f"find {binning_dir} -name 'theta.mrc' > {binning_dir}/theta_files.txt") + + print("\nTo analyze the replicate results with the saved MIP list:") + print(f"python3 scripts/testing/programs/image_replicate_analysis.py --image-list {mip_list_file} --threshold {threshold_value:.3f}") + print("\nOr to extract the threshold value directly from histogram files:") print(f"threshold=$(awk '/Threshold/{{print $NF; exit}}' {hist_filenames[0]})") - print(f"python3 scripts/testing/programs/cistem_test_utils/image_replicate_analysis.py --image-list {mip_list_file} --threshold $threshold") + print(f"python3 scripts/testing/programs/image_replicate_analysis.py --image-list {mip_list_file} --threshold $threshold") # Print the directory where files are saved print(f"\nMIP files saved in: {temp_dir}") From a5a12c051563f007381812a0f4d781ddf4b21b1e Mon Sep 17 00:00:00 2001 From: "copilot-swe-agent[bot]" <198982749+Copilot@users.noreply.github.com> Date: Sat, 24 May 2025 13:22:00 +0000 Subject: [PATCH 18/21] Add --max-threads CLI arg and fix threshold capitalization Co-authored-by: bHimes <3077528+bHimes@users.noreply.github.com> --- scripts/testing/programs/cistem_test_utils/args.py | 4 ++++ .../test_template_reproducibility.py | 14 +++++++++++++- 2 files changed, 17 insertions(+), 1 deletion(-) diff --git a/scripts/testing/programs/cistem_test_utils/args.py b/scripts/testing/programs/cistem_test_utils/args.py index 942e621b0..ee611a6b3 100644 --- a/scripts/testing/programs/cistem_test_utils/args.py +++ b/scripts/testing/programs/cistem_test_utils/args.py @@ -117,6 +117,10 @@ def parse_TM_args(wanted_binary_name): parser.add_argument('--fast-fft', dest='fast_fft', action='store_true', default=True, help='Use FastFFT implementation (default: True)') args_to_check.append('fast_fft') + + parser.add_argument('--max-threads', dest='max_threads', type=int, default=2, + help='Maximum number of threads to use (default: 2)') + args_to_check.append('max_threads') args = parser.parse_args() diff --git a/scripts/testing/programs/match_template/test_template_reproducibility.py b/scripts/testing/programs/match_template/test_template_reproducibility.py index 78e25582a..5357f8e66 100755 --- a/scripts/testing/programs/match_template/test_template_reproducibility.py +++ b/scripts/testing/programs/match_template/test_template_reproducibility.py @@ -62,6 +62,18 @@ def main(): # Parse command-line arguments args = tmArgs.parse_TM_args(wanted_binary_name) + # Print non-default configuration parameters + print("Configuration parameters:") + print(f" max_threads: {args.max_threads}") + if args.fast_fft != True: # Default is True + print(f" fast_fft: {args.fast_fft}") + if hasattr(args, 'gpu_idx') and args.gpu_idx != 0: # Default is 0 + print(f" gpu_idx: {args.gpu_idx}") + if args.cpu: # Default is False + print(" Using CPU implementation") + if args.old_cistem: # Default is False + print(" Using old cisTEM implementation") + # Handle temp directory management options if args.list_temp_dirs: temp_manager.print_temp_dirs() @@ -190,7 +202,7 @@ def main(): print(f"python3 scripts/testing/programs/image_replicate_analysis.py --image-list {mip_list_file} --threshold {threshold_value:.3f}") print("\nOr to extract the threshold value directly from histogram files:") - print(f"threshold=$(awk '/Threshold/{{print $NF; exit}}' {hist_filenames[0]})") + print(f"threshold=$(awk '/threshold/{{print $NF; exit}}' {hist_filenames[0]})") print(f"python3 scripts/testing/programs/image_replicate_analysis.py --image-list {mip_list_file} --threshold $threshold") # Print the directory where files are saved From 2ebdcd46b37f82a2e71387f6df0fd7744b1d43a5 Mon Sep 17 00:00:00 2001 From: "copilot-swe-agent[bot]" <198982749+Copilot@users.noreply.github.com> Date: Sat, 24 May 2025 13:48:41 +0000 Subject: [PATCH 19/21] Implement CLI vs default config comparison and remove problematic file Co-authored-by: bHimes <3077528+bHimes@users.noreply.github.com> --- .../programs/cistem_test_utils/args.py | 20 ++++++++++++++++++- .../match_template/.cistem_temp_dirs.log | 1 - .../test_template_reproducibility.py | 12 ----------- 3 files changed, 19 insertions(+), 14 deletions(-) delete mode 100644 scripts/testing/programs/match_template/.cistem_temp_dirs.log diff --git a/scripts/testing/programs/cistem_test_utils/args.py b/scripts/testing/programs/cistem_test_utils/args.py index ee611a6b3..68641a067 100644 --- a/scripts/testing/programs/cistem_test_utils/args.py +++ b/scripts/testing/programs/cistem_test_utils/args.py @@ -62,7 +62,25 @@ def get_config(args, data_dir: str, ref_number: int, img_number: int): config['result_ignore_n_pixels_from_edge'] = -1 for arg_val in args.args_to_check: - config[arg_val] = getattr(args, arg_val) + # Store the default value for comparison + default_val = config.get(arg_val) + # Get the value from args + arg_val_value = getattr(args, arg_val) + + # Compare with default value if it exists in config and the type matches + if arg_val in config and arg_val_value != default_val: + # Handle string vs numeric comparison + if (isinstance(default_val, (int, float)) and isinstance(arg_val_value, (int, float))) or \ + (isinstance(default_val, str) and isinstance(arg_val_value, str)) or \ + (isinstance(default_val, bool) and isinstance(arg_val_value, bool)): + print(f"User has set {arg_val} to value {arg_val_value}. Changing from default {default_val}") + else: + # Different types or one is None, just print if both are not None + if default_val is not None and arg_val_value is not None: + print(f"User has set {arg_val} to value {arg_val_value}. Changing from default {default_val}") + + # Update the config with the new value + config[arg_val] = arg_val_value config['output_file_prefix'] = os.path.abspath(os.path.join(args.output_file_prefix, config.get('data')[img_number]['img_name'])) os.makedirs(config['output_file_prefix'], exist_ok=True) diff --git a/scripts/testing/programs/match_template/.cistem_temp_dirs.log b/scripts/testing/programs/match_template/.cistem_temp_dirs.log deleted file mode 100644 index 0637a088a..000000000 --- a/scripts/testing/programs/match_template/.cistem_temp_dirs.log +++ /dev/null @@ -1 +0,0 @@ -[] \ No newline at end of file diff --git a/scripts/testing/programs/match_template/test_template_reproducibility.py b/scripts/testing/programs/match_template/test_template_reproducibility.py index 5357f8e66..3cf0781b7 100755 --- a/scripts/testing/programs/match_template/test_template_reproducibility.py +++ b/scripts/testing/programs/match_template/test_template_reproducibility.py @@ -62,18 +62,6 @@ def main(): # Parse command-line arguments args = tmArgs.parse_TM_args(wanted_binary_name) - # Print non-default configuration parameters - print("Configuration parameters:") - print(f" max_threads: {args.max_threads}") - if args.fast_fft != True: # Default is True - print(f" fast_fft: {args.fast_fft}") - if hasattr(args, 'gpu_idx') and args.gpu_idx != 0: # Default is 0 - print(f" gpu_idx: {args.gpu_idx}") - if args.cpu: # Default is False - print(" Using CPU implementation") - if args.old_cistem: # Default is False - print(" Using old cisTEM implementation") - # Handle temp directory management options if args.list_temp_dirs: temp_manager.print_temp_dirs() From 59ae785830eb23fdff8995b14a114edb60333d01 Mon Sep 17 00:00:00 2001 From: "copilot-swe-agent[bot]" <198982749+Copilot@users.noreply.github.com> Date: Sat, 24 May 2025 13:56:16 +0000 Subject: [PATCH 20/21] Simplify CLI argument comparison in args.py Co-authored-by: bHimes <3077528+bHimes@users.noreply.github.com> --- scripts/testing/programs/cistem_test_utils/args.py | 12 +++--------- 1 file changed, 3 insertions(+), 9 deletions(-) diff --git a/scripts/testing/programs/cistem_test_utils/args.py b/scripts/testing/programs/cistem_test_utils/args.py index 68641a067..8b4c62ca5 100644 --- a/scripts/testing/programs/cistem_test_utils/args.py +++ b/scripts/testing/programs/cistem_test_utils/args.py @@ -67,17 +67,11 @@ def get_config(args, data_dir: str, ref_number: int, img_number: int): # Get the value from args arg_val_value = getattr(args, arg_val) - # Compare with default value if it exists in config and the type matches + # Compare with default value if it exists in config and the arg value is not None if arg_val in config and arg_val_value != default_val: - # Handle string vs numeric comparison - if (isinstance(default_val, (int, float)) and isinstance(arg_val_value, (int, float))) or \ - (isinstance(default_val, str) and isinstance(arg_val_value, str)) or \ - (isinstance(default_val, bool) and isinstance(arg_val_value, bool)): + # Only print if both values are not None + if default_val is not None and arg_val_value is not None: print(f"User has set {arg_val} to value {arg_val_value}. Changing from default {default_val}") - else: - # Different types or one is None, just print if both are not None - if default_val is not None and arg_val_value is not None: - print(f"User has set {arg_val} to value {arg_val_value}. Changing from default {default_val}") # Update the config with the new value config[arg_val] = arg_val_value From 0c031615f09446d844631c4e2113f9092b9d6f03 Mon Sep 17 00:00:00 2001 From: "copilot-swe-agent[bot]" <198982749+Copilot@users.noreply.github.com> Date: Sat, 24 May 2025 14:35:52 +0000 Subject: [PATCH 21/21] Remove Python package related files and instructions Co-authored-by: bHimes <3077528+bHimes@users.noreply.github.com> --- scripts/testing/README.md | 26 -------------- .../cistem_test_utils.egg-info/PKG-INFO | 4 --- .../cistem_test_utils.egg-info/SOURCES.txt | 6 ---- .../dependency_links.txt | 1 - .../cistem_test_utils.egg-info/requires.txt | 3 -- .../cistem_test_utils.egg-info/top_level.txt | 1 - scripts/testing/programs/cistem_path_setup.py | 22 ------------ .../cistem_test_utils.egg-info/PKG-INFO | 4 --- .../cistem_test_utils.egg-info/SOURCES.txt | 16 --------- .../dependency_links.txt | 1 - .../cistem_test_utils.egg-info/requires.txt | 3 -- .../cistem_test_utils.egg-info/top_level.txt | 1 - scripts/testing/programs/setup_env.py | 35 ------------------- scripts/testing/setup.py | 14 -------- 14 files changed, 137 deletions(-) delete mode 100644 scripts/testing/cistem_test_utils.egg-info/PKG-INFO delete mode 100644 scripts/testing/cistem_test_utils.egg-info/SOURCES.txt delete mode 100644 scripts/testing/cistem_test_utils.egg-info/dependency_links.txt delete mode 100644 scripts/testing/cistem_test_utils.egg-info/requires.txt delete mode 100644 scripts/testing/cistem_test_utils.egg-info/top_level.txt delete mode 100644 scripts/testing/programs/cistem_path_setup.py delete mode 100644 scripts/testing/programs/cistem_test_utils.egg-info/PKG-INFO delete mode 100644 scripts/testing/programs/cistem_test_utils.egg-info/SOURCES.txt delete mode 100644 scripts/testing/programs/cistem_test_utils.egg-info/dependency_links.txt delete mode 100644 scripts/testing/programs/cistem_test_utils.egg-info/requires.txt delete mode 100644 scripts/testing/programs/cistem_test_utils.egg-info/top_level.txt delete mode 100644 scripts/testing/programs/setup_env.py delete mode 100644 scripts/testing/setup.py diff --git a/scripts/testing/README.md b/scripts/testing/README.md index 34ee336c6..258c202f0 100644 --- a/scripts/testing/README.md +++ b/scripts/testing/README.md @@ -2,32 +2,6 @@ This directory contains test scripts and utilities for testing cisTEM functionality. -## Setup - -To properly use these test scripts, you should install the `cistem_test_utils` package in development mode: - -```bash -# Navigate to the testing directory -cd /path/to/cisTEM/scripts/testing - -# Install the package in development mode -pip install -e . -``` - -## Required Dependencies - -The test scripts require the following Python packages: - -- mrcfile -- toml -- numpy - -These can be installed via: - -```bash -pip install mrcfile toml numpy -``` - ## Running Tests Each test can be run individually from their respective directories. For example: diff --git a/scripts/testing/cistem_test_utils.egg-info/PKG-INFO b/scripts/testing/cistem_test_utils.egg-info/PKG-INFO deleted file mode 100644 index 8694c5e1d..000000000 --- a/scripts/testing/cistem_test_utils.egg-info/PKG-INFO +++ /dev/null @@ -1,4 +0,0 @@ -Metadata-Version: 2.1 -Name: cistem-test-utils -Version: 0.1.0 -Summary: Testing utilities for cisTEM diff --git a/scripts/testing/cistem_test_utils.egg-info/SOURCES.txt b/scripts/testing/cistem_test_utils.egg-info/SOURCES.txt deleted file mode 100644 index 87e6c43d5..000000000 --- a/scripts/testing/cistem_test_utils.egg-info/SOURCES.txt +++ /dev/null @@ -1,6 +0,0 @@ -setup.py -cistem_test_utils.egg-info/PKG-INFO -cistem_test_utils.egg-info/SOURCES.txt -cistem_test_utils.egg-info/dependency_links.txt -cistem_test_utils.egg-info/requires.txt -cistem_test_utils.egg-info/top_level.txt \ No newline at end of file diff --git a/scripts/testing/cistem_test_utils.egg-info/dependency_links.txt b/scripts/testing/cistem_test_utils.egg-info/dependency_links.txt deleted file mode 100644 index 8b1378917..000000000 --- a/scripts/testing/cistem_test_utils.egg-info/dependency_links.txt +++ /dev/null @@ -1 +0,0 @@ - diff --git a/scripts/testing/cistem_test_utils.egg-info/requires.txt b/scripts/testing/cistem_test_utils.egg-info/requires.txt deleted file mode 100644 index cbc2e3e74..000000000 --- a/scripts/testing/cistem_test_utils.egg-info/requires.txt +++ /dev/null @@ -1,3 +0,0 @@ -mrcfile -numpy -toml diff --git a/scripts/testing/cistem_test_utils.egg-info/top_level.txt b/scripts/testing/cistem_test_utils.egg-info/top_level.txt deleted file mode 100644 index 8b1378917..000000000 --- a/scripts/testing/cistem_test_utils.egg-info/top_level.txt +++ /dev/null @@ -1 +0,0 @@ - diff --git a/scripts/testing/programs/cistem_path_setup.py b/scripts/testing/programs/cistem_path_setup.py deleted file mode 100644 index 148e85d6d..000000000 --- a/scripts/testing/programs/cistem_path_setup.py +++ /dev/null @@ -1,22 +0,0 @@ -""" -Simple module to replace annoying_hack.py for cisTEM test scripts. - -This module should be imported at the start of any test script that needs -access to cistem_test_utils. It adds the programs directory to the Python path. - -Usage: - import cistem_path_setup # Add this as the first import in test scripts - - # Then import cistem_test_utils modules as normal - import cistem_test_utils.args as tmArgs -""" - -import os -import sys - -# Get the absolute path to the programs directory (parent of this file) -cistem_programs_path = os.path.dirname(os.path.abspath(__file__)) - -# Add to Python path if not already there -if cistem_programs_path not in sys.path: - sys.path.append(cistem_programs_path) \ No newline at end of file diff --git a/scripts/testing/programs/cistem_test_utils.egg-info/PKG-INFO b/scripts/testing/programs/cistem_test_utils.egg-info/PKG-INFO deleted file mode 100644 index 8694c5e1d..000000000 --- a/scripts/testing/programs/cistem_test_utils.egg-info/PKG-INFO +++ /dev/null @@ -1,4 +0,0 @@ -Metadata-Version: 2.1 -Name: cistem-test-utils -Version: 0.1.0 -Summary: Testing utilities for cisTEM diff --git a/scripts/testing/programs/cistem_test_utils.egg-info/SOURCES.txt b/scripts/testing/programs/cistem_test_utils.egg-info/SOURCES.txt deleted file mode 100644 index bed479471..000000000 --- a/scripts/testing/programs/cistem_test_utils.egg-info/SOURCES.txt +++ /dev/null @@ -1,16 +0,0 @@ -README.md -setup.py -programs/cistem_test_utils/__init__.py -programs/cistem_test_utils/annoying_hack.py -programs/cistem_test_utils/args.py -programs/cistem_test_utils/compare_coordinates.py -programs/cistem_test_utils/make_tmp_runfile.py -programs/cistem_test_utils/re_run_results_on_mip.py -programs/cistem_test_utils/run_job.py -programs/cistem_test_utils/temp_dir_manager.py -programs/cistem_test_utils/threshold_utils.py -programs/cistem_test_utils.egg-info/PKG-INFO -programs/cistem_test_utils.egg-info/SOURCES.txt -programs/cistem_test_utils.egg-info/dependency_links.txt -programs/cistem_test_utils.egg-info/requires.txt -programs/cistem_test_utils.egg-info/top_level.txt \ No newline at end of file diff --git a/scripts/testing/programs/cistem_test_utils.egg-info/dependency_links.txt b/scripts/testing/programs/cistem_test_utils.egg-info/dependency_links.txt deleted file mode 100644 index 8b1378917..000000000 --- a/scripts/testing/programs/cistem_test_utils.egg-info/dependency_links.txt +++ /dev/null @@ -1 +0,0 @@ - diff --git a/scripts/testing/programs/cistem_test_utils.egg-info/requires.txt b/scripts/testing/programs/cistem_test_utils.egg-info/requires.txt deleted file mode 100644 index cbc2e3e74..000000000 --- a/scripts/testing/programs/cistem_test_utils.egg-info/requires.txt +++ /dev/null @@ -1,3 +0,0 @@ -mrcfile -numpy -toml diff --git a/scripts/testing/programs/cistem_test_utils.egg-info/top_level.txt b/scripts/testing/programs/cistem_test_utils.egg-info/top_level.txt deleted file mode 100644 index 7421fb8e4..000000000 --- a/scripts/testing/programs/cistem_test_utils.egg-info/top_level.txt +++ /dev/null @@ -1 +0,0 @@ -cistem_test_utils diff --git a/scripts/testing/programs/setup_env.py b/scripts/testing/programs/setup_env.py deleted file mode 100644 index 3729b0df4..000000000 --- a/scripts/testing/programs/setup_env.py +++ /dev/null @@ -1,35 +0,0 @@ -#!/usr/bin/env python3 -""" -Environmental setup utility for cisTEM test scripts. - -This script modifies the Python sys.path to include the necessary directories -so that test scripts can find and import cistem_test_utils without requiring -annoying_hack.py or special environment variables. -""" - -import os -import sys - -def setup_cistem_env(): - """ - Add the necessary directories to sys.path for cisTEM testing. - - This function is designed to be imported and called at the beginning of - test scripts to ensure all the required modules are available. - """ - # Get the absolute path to the programs directory - script_dir = os.path.dirname(os.path.abspath(__file__)) - - # Check if the directories are already in sys.path - if script_dir not in sys.path: - sys.path.append(script_dir) - - return True - -if __name__ == "__main__": - # If run directly, just report success - if setup_cistem_env(): - print("cisTEM environment setup successfully.") - else: - print("Failed to set up cisTEM environment.") - sys.exit(1) \ No newline at end of file diff --git a/scripts/testing/setup.py b/scripts/testing/setup.py deleted file mode 100644 index 71a80cd82..000000000 --- a/scripts/testing/setup.py +++ /dev/null @@ -1,14 +0,0 @@ -from setuptools import setup, find_packages - -setup( - name="cistem_test_utils", - version="0.1.0", - packages=find_packages(where="programs"), - package_dir={"": "programs"}, - install_requires=[ - "mrcfile", - "toml", - "numpy" - ], - description="Testing utilities for cisTEM", -) \ No newline at end of file