diff --git a/.github/workflows/tests.yml b/.github/workflows/tests.yml index 30ec5f0..0746ed2 100644 --- a/.github/workflows/tests.yml +++ b/.github/workflows/tests.yml @@ -2,9 +2,9 @@ name: Tests on: push: - branches: [ main, dev_1.4, dev_1.7, dev_v0.2.0, dev_v0.2.1] + branches: [ main, dev_1.4, dev_1.7, dev_v0.2.0, dev_v0.2.1, dev_v0.2.2 ] pull_request: - branches: [ main, dev_1.4, dev_1.7, dev_v0.2.0, dev_v0.2.1 ] + branches: [ main, dev_1.4, dev_1.7, dev_v0.2.0, dev_v0.2.1, dev_v0.2.2 ] jobs: tests: @@ -15,7 +15,7 @@ jobs: fail-fast: false matrix: os: [ubuntu-latest, macos-latest, windows-latest] - python-version: ["3.9", "3.10", "3.11", "3.12"] + python-version: ["3.9", "3.10", "3.11", "3.12", "3.13", "3.14"] steps: # 1) Rรฉcupรฉrer le code diff --git a/.gitignore b/.gitignore index 30aa957..a3f52eb 100644 --- a/.gitignore +++ b/.gitignore @@ -14,6 +14,7 @@ __pycache__/ *.pyd *.dll .pytest_tmp/ +tmp/ *venv*/ *.db *.db* @@ -27,6 +28,9 @@ documentation/readthedocs/_build/ # PACKAGE DATA # ------------------------------- +# Dev TODOs +future_todos.md + # Ignore everything inside OUTPUT, but keep the folder itself src/shinier/data/OUTPUT/* !src/shinier/data/OUTPUT/.gitkeep diff --git a/README.md b/README.md index 50bce5c..8428608 100644 --- a/README.md +++ b/README.md @@ -1,4 +1,4 @@ -# ๐ŸŒŸ SHINIER +# SHINIER ```text โ–ˆโ–ˆโ–ˆโ–ˆโ–ˆโ–ˆโ–ˆโ•—โ–ˆโ–ˆโ•— โ–ˆโ–ˆโ•—โ–ˆโ–ˆโ•—โ–ˆโ–ˆโ–ˆโ•— โ–ˆโ–ˆโ•—โ–ˆโ–ˆโ•—โ–ˆโ–ˆโ–ˆโ–ˆโ–ˆโ–ˆโ–ˆโ•—โ–ˆโ–ˆโ–ˆโ–ˆโ–ˆโ–ˆโ•— โ–ˆโ–ˆโ•”โ•โ•โ•โ•โ•โ–ˆโ–ˆโ•‘ โ–ˆโ–ˆโ•‘โ–ˆโ–ˆโ•‘โ–ˆโ–ˆโ–ˆโ–ˆโ•— โ–ˆโ–ˆโ•‘โ–ˆโ–ˆโ•‘โ–ˆโ–ˆโ•”โ•โ•โ•โ•โ•โ–ˆโ–ˆโ•”โ•โ•โ–ˆโ–ˆโ•— @@ -15,26 +15,26 @@ [![Tests](https://github.com/Charestlab/shinier/actions/workflows/tests.yml/badge.svg)](https://github.com/Charestlab/shinier/actions/workflows/tests.yml) --- -## ๐ŸŽฏ Overview +## Overview SHINIER is a modern Python implementation of SHINE (Spectrum, Histogram, and Intensity Normalization and Equalization), originally developed in MATLAB by Willenbockel etย al., 2010. It provides precise control over luminance, contrast, histograms, and spectral content across large image sets for well-calibrated visual experiments. ### Key Features and Improvements -- ๐ŸŽจ **Color Processing** โ€” New modes for color image control with modern color-space standards (Rec.601 / Rec.709 / Rec.2020). -- ๐Ÿ–ผ๏ธ **Dithering Support** โ€” Reduces quantization artifacts and enhances output image quality. -- โšก **Optimized Performance** โ€” Efficient memory management and faster processing for large image sets (optional Cython/C++ convolution core). -- ๐Ÿ•ฐ **Legacy Mode** โ€” Ensures full backward compatibility with MATLABโ€™s original SHINE toolbox. -- ๐Ÿ”ข **High-Precision Arithmetic** โ€” Computations in floating-point precision rather than 8-bit integer space, minimizing rounding errors in multi-stage processing. -- ๐Ÿ“ฆ **Object-Oriented Design** โ€” Modular, extensible architecture with a clean Python API. -- ๐Ÿ˜€ **User-Friendly CLI** โ€” Guided, prompt-based interface for users who prefer not to write code. +- **Color Processing** โ€” New modes for color image control with modern color-space standards (Rec.601 / Rec.709 / Rec.2020). +- **Dithering Support** โ€” Reduces quantization artifacts and enhances output image quality. +- **Optimized Performance** โ€” Efficient memory management and faster processing for large image sets (optional Cython/C++ convolution core). +- **Legacy Mode** โ€” Ensures full backward compatibility with MATLABโ€™s original SHINE toolbox. +- **High-Precision Arithmetic** โ€” Computations in floating-point precision rather than 8-bit integer space, minimizing rounding errors in multi-stage processing. +- **Object-Oriented Design** โ€” Modular, extensible architecture with a clean Python API. +- **User-Friendly CLI** โ€” Guided, prompt-based interface for users who prefer not to write code. For detailed technical documentation (algorithms, numerical choices, and MATLAB vs Python behavior), see [`documentation/documentation.md`](documentation/documentation.md). --- -## ๐Ÿš€ Quick Start +## Quick Start ### Installation @@ -68,7 +68,7 @@ print("shinier version:", getattr(shinier, "__version__", "unknown")) ``` -### ๐Ÿ˜€ **User-friendly Interface** +### User-friendly Interface Call the following bash command to quickly start using the interactive CLI. ```bash shinier --show_results --image_index=1 @@ -77,7 +77,7 @@ shinier --show_results --image_index=1 CLI demo

-### ๐Ÿงฉ Example in Python +### Example in Python Run the following python code to make sure the package is running properly. ```python from shinier import Options, ImageDataset, ImageProcessor, utils @@ -102,7 +102,7 @@ Change the mode number (e.g. `opt = Options(mode=3)`) to change image processing | 6 | `hist_match โ†’ spec_match` | Histogram, then spectrum | | 7 | `sf_match โ†’ hist_match` | Spatial frequency, then histogram | | 8 | `spec_match โ†’ hist_match` (default) | Spectrum, then histogram (recommended) | -| 9 | `dithering` | Dithering only | +| 9 | `ie_methods` or `dithering` | Standalone per-image transform (histogram-derived enhancement or dithering) | Below is an example of results obtained using mode 5 with joint histogram equalization and spatial frequency normalization.

@@ -110,7 +110,7 @@ Below is an example of results obtained using mode 5 with joint histogram equali

--- -## ๐Ÿ›๏ธ **Technical information** +## **Technical information** Documentation: 1. [Package Overview](documentation/documentation.md#overview) @@ -128,29 +128,30 @@ Documentation: 13. [Additional Resources](documentation/documentation.md#additional-resources) --- -## ๐Ÿ“š Citing +## Citing If you use **SHINIER**, please cite both of these articles: ### References -- Salvas-Hรฉbert, M., Dupuis-Roy, N., Landry, C., Charest, I., & Gosselin, F. (2026). *SHINIER: An Open-Source Python Package for Controlling Low-Level Image Properties* -- Willenbockel, V., Sadr, J., Fiset, D., Horne, G. O., Gosselin, F., & Tanaka, J. W. (2010). Controlling low-level image properties: The SHINE toolbox. *Behavior Research Methods, 42*(3), 671โ€“684. https://doi.org/10.3758/BRM.42.3.671 +- Salvas-Hรฉbert, M., Dupuis-Roy, N., Landry, C., Charest, I., & Gosselin, F. (2026). SHINIER: An open-source Python package for controlling low-level image properties. *SoftwareX*, *35*, Article 102884. https://doi.org/10.1016/j.softx.2026.102884 +- Willenbockel, V., Sadr, J., Fiset, D., Horne, G. O., Gosselin, F., & Tanaka, J. W. (2010). Controlling low-level image properties: The SHINE toolbox. *Behavior Research Methods*, *42*(3), 671โ€“684. https://doi.org/10.3758/BRM.42.3.671 --- -## ๐Ÿค Contributing +## Contributing See [CONTRIBUTING.md](documentation/contributing.md) for guidelines (coding standards, tests, docs, and PR flow). --- -## ๐Ÿ“„ License +## License See [LICENSE](LICENSE) for more information. --- -## ๐Ÿ› ๏ธ Troubleshooting +## Troubleshooting - No compiler available: install a C/C++ toolchain or proceed with the NumPy fallback (slower). +- Compiled extension fails to load (e.g., NumPy too old / ABI mismatch): SHINIER warns and uses the NumPy fallback; upgrade NumPy (and pip/setuptools/wheel) and reinstall to restore the compiled core. - Import errors after upgrade: try pip install --upgrade pip setuptools wheel and reinstall. - Windows build issues: ensure MSVC Build Tools are installed and on PATH. ---

Code developed by Nicolas Dupuis-Roy and Mathias Salvas-Hรฉbert
- Version 0.2.0 - Complete technical documentation + Version 0.2.2 - Complete technical documentation

diff --git a/documentation/contributing.md b/documentation/contributing.md index a89a82a..4589e35 100644 --- a/documentation/contributing.md +++ b/documentation/contributing.md @@ -40,7 +40,7 @@ By participating, you agree to uphold a standard of professional, inclusive, and ## Development Setup -> **Python:** >=3.9, <3.13 +> **Python:** >=3.9, <3.15 > **OS:** macOS / Linux / Windows > **Optional:** C/C++ toolchain for the Cython-compiled `_cconvolve` extension (speeds up convolution) @@ -238,6 +238,6 @@ pytest -m unit_tests # Validation Tests (VERY LONG AND SLOW) and markers pytest -q -m validation_tests ``` -#### ๐Ÿ›๏ธ Test README +#### Test README For complete testing procedures (markers, sharding, replay/debug workflow), see the dedicated Testing Guide in `tests/README.md`. diff --git a/documentation/demos.md b/documentation/demos.md index 3c63bfb..a3fb614 100644 --- a/documentation/demos.md +++ b/documentation/demos.md @@ -140,7 +140,7 @@ MASKS_FOLDER = "path" --- -### Mode 1 โ€“ `lum_match` +#### Mode 1 โ€“ `lum_match` ```python """ @@ -168,7 +168,7 @@ opts = Options( --- -### Mode 2 โ€“ `hist_match` +#### Mode 2 โ€“ `hist_match` ```python """ @@ -193,7 +193,7 @@ opts = Options( --- -### Mode 3 โ€“ `sf_match` +#### Mode 3 โ€“ `sf_match` ```python """ @@ -216,7 +216,7 @@ opts = Options( --- -### Mode 4 โ€“ `spec_match` +#### Mode 4 โ€“ `spec_match` ```python """ @@ -240,7 +240,7 @@ opts = Options( --- -### Mode 5 โ€“ `hist_match` โ†’ `sf_match` +#### Mode 5 โ€“ `hist_match` โ†’ `sf_match` ```python """ @@ -249,7 +249,7 @@ Mode 5 (hist_match & sf_match): histogram matching followed by rotational Example use case: Histogram specification with noise is applied (legacy method), then rotational Fourier spectra are aligned. No rescaling is performed - afterwards,to preserve the luminance distribution imposed by histogram + afterwards, to preserve the luminance distribution imposed by histogram matching. """ opts = Options( @@ -264,7 +264,7 @@ opts = Options( --- -### Mode 6 โ€“ `hist_match` โ†’ `spec_match` +#### Mode 6 โ€“ `hist_match` โ†’ `spec_match` ```python """ @@ -285,7 +285,7 @@ opts = Options( --- -### Mode 7 โ€“ `sf_match` โ†’ `hist_match` +#### Mode 7 โ€“ `sf_match` โ†’ `hist_match` ```python """ @@ -308,7 +308,7 @@ opts = Options( --- -### Mode 8 โ€“ `spec_match` โ†’ `hist_match` +#### Mode 8 โ€“ `spec_match` โ†’ `hist_match` ```python """ @@ -332,11 +332,11 @@ opts = Options( --- -### Mode 9 โ€“ Standalone Per-Image Transform +#### Mode 9 โ€“ Standalone Per-Image Transform Mode 9 applies a standalone transform to each image independently โ€” no inter-image target is computed. The operation is selected via `standalone_op`. -#### 9a โ€“ Image Enhancement (`standalone_op="ie_methods"`) +**9a โ€” Image enhancement** (`standalone_op="ie_methods"`) ```python """ @@ -355,7 +355,7 @@ opts = Options( ) ``` -#### 9b โ€“ Dithering only (`standalone_op="dithering"`) +**9b โ€” Dithering only** (`standalone_op="dithering"`) ```python """ @@ -373,7 +373,7 @@ opts = Options( --- -### Example 10 โ€“ Mode 2 + extra parameters +#### Example 10 โ€“ Mode 2 + extra parameters ```python """ @@ -381,7 +381,7 @@ Example 10 (mode 2 + non-mode-specific parameters): to show the other parameters Example use case: hist_matching using Coltuc, Bolon & Chassery (2006) exact histogram specification. Target histogram will be the average from all the - images (default), no SSIM optimization (Avanki, 2009). + images (default), no SSIM optimization (Avanaki, 2009). The masks loaded from masks_folder are used for figure-ground separation (whole_image = 3). The images will be transformed to grayscale (1 channel), @@ -405,7 +405,7 @@ opts = Options( --- -### Example 11 โ€“ Preserve colors with `constrain_image_chrominance` +#### Example 11 โ€“ Preserve colors with `constrain_image_chrominance` ```python """ diff --git a/documentation/documentation.md b/documentation/documentation.md index 60ad522..ba9ab3c 100644 --- a/documentation/documentation.md +++ b/documentation/documentation.md @@ -15,26 +15,31 @@ # Documentation -## ๐Ÿ“‹ Table of Contents +## Table of Contents 1. [Overview](#overview) 2. [Package Architecture](#package-architecture) 3. [MATLAB vs Python Differences](#matlab-vs-python-differences) 4. [Detailed Processing Modes](#detailed-processing-modes) -5. [Main Classes](#main-classes) -6. [Visualization Functions](#visualization-functions) -7. [Implemented Algorithms](#implemented-algorithms) -8. [Memory Management and Performance](#memory-management-and-performance) -9. [Testing and Validation](#testing-and-validation) +5. [Border Artifacts and FFT Padding](#border-artifacts-and-fft-padding) +6. [Main Classes](#main-classes) +7. [Visualization Functions](#visualization-functions) +8. [Implemented Algorithms](#implemented-algorithms) +9. [Memory Management and Performance](#memory-management-and-performance) +10. [Testing and Validation](#testing-and-validation) +11. [Troubleshooting and Optimization](#troubleshooting-and-optimization) +12. [Additional Resources](#additional-resources) --- -## ๐ŸŽฏ Overview +## Overview **SHINIER** is a modern Python implementation of the **SHINE** (Spectrum, Histogram, and Intensity Normalization, Equalization, and Refinements) toolbox, originally developed in MATLAB by [Willenbockel et al. (2010)](https://doi.org/10.3758/BRM.42.3.671). This new version implemented new options (e.g., color management, dithering and [Coltuc, Bolon & Chassery (2006)](https://www.cin.ufpe.br/~if751/projetos/artigos/Exact%20Histogram%20Specification.pdf) exact histogram specification algorithm) and refined the previous ones. -**Reference** : [Willenbockel, V., Sadr, J., Fiset, D., Horne, G. O., Gosselin, F., & Tanaka, J. W. (2010). Controlling low-level image properties: The SHINE toolbox. *Behavior Research Methods*, 42(3), 671-684.](https://doi.org/10.3758/BRM.42.3.671) +**References**: +- [Salvas-Hรฉbert, M., Dupuis-Roy, N., Landry, C., Charest, I., & Gosselin, F. (2026). SHINIER: An open-source Python package for controlling low-level image properties. *SoftwareX*, *35*, Article 102884.](https://doi.org/10.1016/j.softx.2026.102884) +- [Willenbockel, V., Sadr, J., Fiset, D., Horne, G. O., Gosselin, F., & Tanaka, J. W. (2010). Controlling low-level image properties: The SHINE toolbox. *Behavior Research Methods*, *42*(3), 671โ€“684.](https://doi.org/10.3758/BRM.42.3.671) ### Main Objectives - **Compatibility**: Maintain compatibility with the original MATLAB implementation @@ -45,7 +50,7 @@ --- -## ๐Ÿ—๏ธ Package Architecture +## Package Architecture ### Module Structure @@ -73,7 +78,7 @@ shinier/src --- -## ๐Ÿ”ฌ MATLAB vs Python Differences +## MATLAB vs Python Differences ### 1. **Rounding Operators** @@ -258,9 +263,11 @@ At half-integer boundaries this shifts ~0.1 % of pixels by ยฑ1 gray level. There --- -## โš™๏ธ Detailed Processing Modes +## Detailed Processing Modes -### Mode 1: Luminance Matching Only +### **Pixel-based matching (Modes 1โ€“2)** + +#### Mode 1: Luminance Matching Only ```python mode = 1 # lum_match only ``` @@ -278,10 +285,15 @@ mode = 1 # lum_match only - `target_lum`: Tuple `(mean, std)` where `mean โˆˆ [0, 255]` or `None`, and `std โˆˆ [0, +โˆž)` or `None`. `0` uses the dataset average for that statistic, so `(0, 20)` uses the average mean and a std of 20, while `(100, 0)` uses a mean of 100 and the average std. - `safe_lum_match`: If True, automatically adjusts `(target_mean, target_std)` to keep all pixel values within [0, 255] (values may differ slightly from the requested target) -### Mode 2: Histogram Matching Only +#### Mode 2: Histogram Matching Only ```python mode = 2 # hist_match only ``` + +**What it does (sliding-puzzle analogy).** A histogram is only an *inventory* of shades โ€” how many pixels are black, dark gray, gray, or white โ€” with no information about *where* each one goes. The same inventory can be arranged into completely different pictures: below, the scrambled panel and both faces share the **exact same histogram** (13 / 9 / 106 / 96 pixels). Like a sliding puzzle, histogram matching redistributes an image's pixels to reproduce a *target's* inventory of shades, while preserving each pixel's brightness rank so the original structure is kept. + +![](figures/sliding_puzzle.png) + **Available Algorithms:** - **Exact specification** (`hist_specification=0`): [Coltuc, Bolon & Chassery (2006)]((https://www.cin.ufpe.br/~if751/projetos/artigos/Exact%20Histogram%20Specification.pdf)) algorithm - **Specification with noise** (`hist_specification=1`): Legacy version with noise addition @@ -291,7 +303,9 @@ mode = 2 # hist_match only - `hist_iterations`: Number of iterations (default: 10) - `step_size`: Step size (default: 34) -### Mode 3: Spatial Frequency Matching Only +### **Spatial-frequency-based matching (Modes 3โ€“4)** + +#### Mode 3: Spatial Frequency Matching Only ```python mode = 3 # sf_match only ``` @@ -307,7 +321,7 @@ Equalizes the **mean amplitude per spatial frequency** across images โ€” i.e., t --- -### Mode 4: Spectrum Matching Only +#### Mode 4: Spectrum Matching Only ```python mode = 4 # spec_match only ``` @@ -321,7 +335,7 @@ Equalizes the **amplitude at every spatial frequency and orientation** across im 4. Apply rescaling per `options.rescaling` and clip to the valid intensity range. 5. Store float255 outputs and optionally log or visualize spectrum-related diagnostics. -### Composite Modes (5-8) +### **Composite modes (Modes 5โ€“8)** #### Mode 5: Histogram + Spatial Frequency ```python @@ -347,7 +361,9 @@ mode = 8 # spec_match โ†’ hist_match - Composite modes use temporary floating-point precision - Reduces rounding errors in multi-step calculations -### Mode 9: Standalone Per-Image Transform +### **Standalone transforms (Mode 9)** + +#### Mode 9: Standalone per-image transform ```python mode = 9 # standalone_op = "ie_methods" or "dithering" @@ -360,16 +376,15 @@ dithering = 1 # Dithering method (0 = none, 1 = Noisy-bit, 2 = Floyd-Steinberg) Applies a standalone transform to each image independently โ€” no inter-image target is computed. ---- - -### Histogram Equalization: Exact Specification and Histogram-Derived Remapping +#### Histogram equalization: exact specification and histogram-derived remapping Histogram equalization can be achieved through **Exact Histogram Specification (EHS)** using a flat, uniform target histogram (`target_hist="equal"`, `mode=2` or modes 5โ€“8). Pixels are individually ranked and assigned to target bins, allowing the output to exactly match the feasible discrete uniform histogram. SHINIER also provides **histogram-derived methods** (`mode=9`, `standalone_op="ie_methods"`), including `classic_he`, `tidhe`, and `rdfhe`. These methods compute gray-level mappings from the image histogram or CDF. Because identical input intensities receive the same output value, the resulting histogram is generally only approximately uniform. + --- -### Border Artifacts and FFT Padding +## Border Artifacts and FFT Padding **Why border artifacts occur** The Fourier transform implicitly treats an image as if it repeats infinitely in all directions: the left edge connects to the right, and the top to the bottom. When opposite edges differ, this creates artificial discontinuities that introduce unwanted high-frequency energy (*spectral leakage*), sometimes visible as ringing or edge artifacts after reconstruction. @@ -397,7 +412,7 @@ Images can optionally be padded before computing the FFT, then cropped back to t --- -## ๐Ÿ›๏ธ Main Classes +## Main Classes ### `Converter` Encapsulates color-space conversions and transfer functions for Rec.601/709/2020. ```python @@ -563,14 +578,6 @@ class ImageProcessor: ``` -## Additional Resources - -- For a detailed description of the available options, see the `Options` class in `Options.py`; each parameter lists its purpose, allowed values, and default. -- For algorithmic details and a walkthrough of processing steps, see the `ImageProcessor` class in `ImageProcessor.py`. -- For color management and gamut-control strategies, see the `GamutControl` class in `color/GamutControl.py`. Interactive visual examples are available at [shinier-web examples](https://charestlab.github.io/shinier-web/). - ---- - ## Visualization Functions These helpers are implemented in `src/shinier/utils.py`. @@ -690,12 +697,14 @@ masked_images = masker.apply_all(stim_arr) mask_from_gui = masker.interactive_mask(image) ``` +![Dynamic StimulusMasker GUI demo](readthedocs/_static/dynamic_stim_masker.gif) + --- -## ๐Ÿงฎ Implemented Algorithms +## Implemented Algorithms -### 1. Exact Histogram Specification +### 1. **Exact Histogram Specification** **Reference:** [Coltuc, D., Bolon, P., & Chassery, J. M. (2006). Exact histogram specification. *IEEE Transactions on Image Processing*, 15(5), 1143-1152.](https://www.cin.ufpe.br/~if751/projetos/artigos/Exact%20Histogram%20Specification.pdf) **Algorithm:** @@ -704,7 +713,7 @@ mask_from_gui = masker.interactive_mask(image) 3. Create mapping table based on CDFs 4. Apply mapping pixel by pixel -### 2. SSIM Optimization for Histogram +### 2. **SSIM Optimization for Histogram** **Reference:** [Avanaki, A. N. (2009). Exact histogram specification for digital images using a variational approach. *Journal of Visual Communication and Image Representation*, 20(7), 505-515.](https://link.springer.com/article/10.1007/s10043-009-0119-z) **Algorithm:** @@ -712,7 +721,7 @@ mask_from_gui = masker.interactive_mask(image) 2. Successive iterations with SSIM-based adjustment 3. Step size optimization for fast convergence -### 3. Floyd-Steinberg Dithering +### 3. **Floyd-Steinberg Dithering** **Reference:** Floyd, R. W., & Steinberg, L. (1976). An adaptive algorithm for spatial grey scale. **Algorithm:** @@ -720,7 +729,7 @@ mask_from_gui = masker.interactive_mask(image) 2. Calculate quantization error for each pixel 3. Distribute error to neighboring pixels with different weights. -### 4. Noisy Bit Dithering +### 4. **Noisy Bit Dithering** **Reference:** [Allard, R., & Faubert, J. (2008). The noisy-bit method for digital halftoning. *Journal of the Optical Society of America A*, 25(8), 1980-1989.](https://link.springer.com/article/10.3758/BRM.40.3.735) **Algorithm:** @@ -728,7 +737,7 @@ mask_from_gui = masker.interactive_mask(image) 2. Quantize with rounding 3. Preserve overall image structure -### 5. Classic Global Histogram Equalization (Classic HE) +### 5. **Classic Global Histogram Equalization (Classic HE)** **Algorithm:** 1. Compute the intensity histogram of the image @@ -738,7 +747,7 @@ mask_from_gui = masker.interactive_mask(image) The output histogram is approximately flat over [0, 255]. Maximizes contrast globally but can over-enhance noise on natural images. Implemented by `shinier.utils.classic_he_gray`. -### 6. Tripartite Image Decomposition-Based Histogram Equalization (TIDHE) +### 6. **Tripartite Image Decomposition-Based Histogram Equalization (TIDHE)** **Reference:** [Rahman, H., & Shimamura, T. (2026). Tripartite image decomposition-based histogram equalization to enhance slightly low-contrast and low-contrast images. ICIC Express Letters, 20(3), 321-332.](https://doi.org/10.24507/icicel.20.03.321) **Algorithm:** @@ -747,7 +756,7 @@ The output histogram is approximately flat over [0, 255]. Maximizes contrast glo 3. Clip each sub-histogram at the average of its mean and median to control enhancement rate (Eqs. 3-5) 4. Equalize each sub-band independently via its clipped CDF (Eqs. 9-11) -### 7. Recursive Dualistic Fuzzy Histogram Equalization (RDFHE) +### 7. **Recursive Dualistic Fuzzy Histogram Equalization (RDFHE)** **Reference:** [Rahman, H., Mostofa, S., Akter, T., & Rashedunnabi, A. H. M. (2026, April). Efficient enhancement of images using recursive dualistic fuzzy histogram equalization. In *2026 IEEE 2nd International Conference on Quantum Photonics, Artificial Intelligence & Networking (QPAIN)* (pp. 1โ€“6). IEEE.](https://doi.org/10.1109/QPAIN69676.2026.11546014) **Algorithm:** @@ -758,7 +767,7 @@ The output histogram is approximately flat over [0, 255]. Maximizes contrast glo Implemented by `shinier.utils.rdfhe_gray`. -### 8. Nonlinear Fuzzificationโ€“Linear Defuzzification-Based ICE (NFLDICE) +### 8. **Nonlinear Fuzzificationโ€“Linear Defuzzification-Based ICE (NFLDICE)** **Reference:** [Rahman, H. (2025). A Time-Efficient and Effective Image Contrast Enhancement Technique Using Fuzzification and Defuzzification. In *Proceedings of Trends in Electronics and Health Informatics* (Lecture Notes in Networks and Systems, vol. 1034, pp. 45โ€“58). Springer.](https://doi.org/10.1007/978-981-97-3937-0_4) Unlike the fuzzy-histogram methods (DFHE/RDFHE), NFLDICE is a fuzzy **set-theoretic** technique: it operates directly on gray levels as fuzzy sets rather than on the histogram. @@ -770,7 +779,7 @@ Unlike the fuzzy-histogram methods (DFHE/RDFHE), NFLDICE is a fuzzy **set-theore Reference parameters: `B=10`, `E_l=5`, `P_l=127.5`, `L=256`. Implemented by `shinier.utils.nfldice_gray`. -### 9. Bi-Entropy Curve Equalization (BETCE) +### 9. **Bi-Entropy Curve Equalization (BETCE)** **Reference:** [Rahman, H. (2025). Bi-Entropy Curve Equalization for Enhancement of Images. In *2025 7th International Conference on Electrical Information and Communication Technology (EICT)* (pp. 1โ€“6). IEEE.](https://doi.org/10.1109/EICT68394.2025.11355632) BETCE is a state-of-the-art curve-based algorithm for very low-contrast grayscale images. It replaces the image histogram with an entropy curve, partitions that curve into lower and upper sub-curves, and equalizes each sub-curve independently. @@ -783,7 +792,7 @@ BETCE is a state-of-the-art curve-based algorithm for very low-contrast grayscal Implemented by `shinier.utils.betce_gray`. -### 10. Sakaguchi-Type Function-Based Cost-Effective Filtering (SFCEF) +### 10. **Sakaguchi-Type Function-Based Cost-Effective Filtering (SFCEF)** **Reference:** [Rahman, H., Sugiura, Y., & Shimamura, T. (2025). Enhancement of low-light images using Sakaguchi-type function-based cost-effective filtering. *Pattern Analysis and Applications*, 28, 193.](https://doi.org/10.1007/s10044-025-01578-8) SFCEF is a state-of-the-art filtering-based algorithm for low-light grayscale images. It builds one 3x3 convolution filter from coefficient bounds of a Sakaguchi/Gegenbauer geometric function class, then filters the input image directly. @@ -800,7 +809,7 @@ Implemented by `shinier.utils.sfcef_gray`. --- -## ๐Ÿ’พ Memory Management and Performance +## Memory Management and Performance ### Memory Conservation Mode (`conserve_memory=True`) @@ -829,7 +838,7 @@ class ImageListIO: --- -## ๐Ÿงช Testing and Validation +## Testing and Validation ### Unit Tests @@ -856,19 +865,22 @@ See `tests/README.md` for a detailed breakdown. `ImageEnhancement_validation_test.py` validates each image-enhancement algorithm (TIDHE, RDFHE, NFLDICE, BETCE, SFCEF) against pixel-exact MATLAB reference outputs stored as SHA-256 hashes in `tests/assets/image_enhancement_matlab_sha256.json`. SFCEF uses a pixel-difference bound (`max_diff โ‰ค 1`) instead of exact hash equality due to FMA-induced rounding differences between MATLAB and NumPy. ---- +### MATLAB SHINE Comparison -## ๐Ÿ“š Usage Examples +A standalone tool benchmarks SHINIER directly against the original [MATLAB SHINE toolbox](http://www.mapageweb.umontreal.ca/gosselif/SHINE/) across processing modes 1โ€“8. +It compares three implementations โ€” `matlab_shine` (the original toolbox), `shinier_legacy` (`legacy_mode=True`, MATLAB-compatible behavior), and `shinier_modern_gray` (SHINIER defaults on grayscale) โ€” in two stages: pixel differences between saved outputs, and distances to shared fixed targets (histogram and spectrum), each measured in the implementation's own processing domain. -- [See](documentation/demos.ipynb) `demos.ipynb` [in the documentation folder for](documentation/demos.ipynb): - - Coding usage - - Interactive CLI usage +```bash +# Requires MATLAB and the SHINE toolbox +bash tests/tools/run_matlab_shine_comparison.sh +``` -Examples in here have been intentionally minimized; please open the demos for more examples. +Results are written as CSV files under `tmp/matlab_shine_comparison/` and summarized in terminal tables. See `tests/README.md` and the docstring of `tests/tools/matlab_shine_comparison.py` for details. --- -## ๐Ÿ”ง Troubleshooting and Optimization + +## Troubleshooting and Optimization ### Common Issues @@ -914,7 +926,7 @@ options = Options( ``` **Scientific Rationale:** -Composite modes (5-8) apply **two sequential transformations** (e.g., spectrum matching followed by histogram matching). Because each transformation modifies the image in ways that can partially undo the effects of the other, a **single pass rarely yields convergence**. As detailed in the original [SHINE documentation](documentation/papers/Controlling%20low-level%20image%20properties:%20The%20SHINE%20toolbox.pdf), **iterative application** of both steps allows the algorithm to progressively minimize residual discrepancies between the desired luminance distribution and spectral amplitude structure. +Composite modes (5-8) apply **two sequential transformations** (e.g., spectrum matching followed by histogram matching). Because each transformation modifies the image in ways that can partially undo the effects of the other, a **single pass rarely yields convergence**. As detailed in the original [SHINE documentation](../_static/shine_toolbox.pdf), **iterative application** of both steps allows the algorithm to progressively minimize residual discrepancies between the desired luminance distribution and spectral amplitude structure. 1. **Sequential Processing**: Each cycle compensates for the distortions introduced by the preceding transformation (e.g., histogram adjustment altering spectral power). 2. **Convergence**: Repeated alternation drives both properties toward their joint target values. @@ -922,7 +934,23 @@ Composite modes (5-8) apply **two sequential transformations** (e.g., spectrum m --- + +## Additional Resources + +The examples in this documentation are intentionally minimized. For more **complete usage examples**, see `demos.ipynb` in the documentation folder: + +- Coding usage +- Interactive CLI usage + +For a **detailed description** of the available **options**, see the `Options` class in `Options.py`; each parameter lists its purpose, allowed values, and default. + +For **algorithmic details** and a walkthrough of processing steps, **see** the `ImageProcessor` class in `ImageProcessor.py`. + +For **color management** and **gamut-control strategies**, see the `GamutControl` class in `color/GamutControl.py`. Interactive **visual examples** are available at [shinier-web examples](https://charestlab.github.io/shinier-web/). + +--- +

Code developed by Nicolas Dupuis-Roy and Mathias Salvas-Hรฉbert
- Version 0.2.0 - Complete technical documentation + Version 0.2.2 - Complete technical documentation

diff --git a/documentation/figures/sliding_puzzle.png b/documentation/figures/sliding_puzzle.png new file mode 100644 index 0000000..bfa46b5 Binary files /dev/null and b/documentation/figures/sliding_puzzle.png differ diff --git a/documentation/readthedocs/_static/dynamic_stim_masker.gif b/documentation/readthedocs/_static/dynamic_stim_masker.gif new file mode 100644 index 0000000..c7a9b87 Binary files /dev/null and b/documentation/readthedocs/_static/dynamic_stim_masker.gif differ diff --git a/documentation/readthedocs/_static/shine_toolbox.pdf b/documentation/readthedocs/_static/shine_toolbox.pdf new file mode 100644 index 0000000..9a1c299 Binary files /dev/null and b/documentation/readthedocs/_static/shine_toolbox.pdf differ diff --git a/documentation/readthedocs/_static/shinier.pdf b/documentation/readthedocs/_static/shinier.pdf new file mode 100644 index 0000000..b02eccb Binary files /dev/null and b/documentation/readthedocs/_static/shinier.pdf differ diff --git a/documentation/readthedocs/conf.py b/documentation/readthedocs/conf.py index f91a09c..e929630 100644 --- a/documentation/readthedocs/conf.py +++ b/documentation/readthedocs/conf.py @@ -16,7 +16,7 @@ try: from shinier import __version__ except Exception: - __version__ = "0.2.0" + __version__ = "0.2.2" version = __version__ release = __version__ diff --git a/documentation/readthedocs/index.md b/documentation/readthedocs/index.md index 1199210..fee96d5 100644 --- a/documentation/readthedocs/index.md +++ b/documentation/readthedocs/index.md @@ -37,6 +37,7 @@ contributors_guide/testing :hidden: Project Links +SHINIER PDF License ``` diff --git a/documentation/readthedocs/paper.md b/documentation/readthedocs/paper.md new file mode 100644 index 0000000..d544142 --- /dev/null +++ b/documentation/readthedocs/paper.md @@ -0,0 +1,20 @@ +# SHINIER PDF + +The SHINIER PDF can be opened directly here: + +[Open or download the SHINIER PDF](_static/shinier.pdf) + +```{raw} html + +``` diff --git a/documentation/readthedocs/project-links.md b/documentation/readthedocs/project-links.md index 01740bb..4a61100 100644 --- a/documentation/readthedocs/project-links.md +++ b/documentation/readthedocs/project-links.md @@ -2,5 +2,6 @@ Useful external links for SHINIER: +- Article: [SHINIER (ScienceDirect)](https://www.sciencedirect.com/science/article/pii/S2352711026003754) - PyPI: - GitHub: diff --git a/pyproject.toml b/pyproject.toml index 5cd2f0e..b3b9040 100644 --- a/pyproject.toml +++ b/pyproject.toml @@ -1,10 +1,10 @@ [build-system] -requires = ["setuptools>=77", "wheel", "numpy>=1.22"] +requires = ["setuptools>=77", "wheel", "numpy>=2.0"] build-backend = "setuptools.build_meta" [project] name = "shinier" -version = "0.2.1" +version = "0.2.2" description = "Python port of the SHINE toolbox with added options (color management, dithering, EHS), optimized for large image sets." readme = "README.md" license = "BSD-3-Clause" @@ -13,12 +13,11 @@ authors = [ { name = "Nicolas Dupuis-Roy", email = "nicolas.dupuis.roy@umontreal.ca" }, { name = "Mathias Salvas-Hรฉbert", email = "mathias.salvas-hebert@umontreal.ca" } ] -requires-python = ">=3.9,<3.13" +requires-python = ">=3.9,<3.15" dependencies = [ - "numpy>=1.22,<2.1", + "numpy>=1.23", "Pillow>=9.0.0", "tqdm>=4.60", - "cython>=3.0", "matplotlib>=3.9.0", "pydantic>=2.12", ] @@ -27,6 +26,8 @@ classifiers = [ "Programming Language :: Python :: 3.10", "Programming Language :: Python :: 3.11", "Programming Language :: Python :: 3.12", + "Programming Language :: Python :: 3.13", + "Programming Language :: Python :: 3.14", "Operating System :: OS Independent", ] @@ -38,6 +39,7 @@ dev = [ "pytest>=8.3.0", "wheel>=0.44.0", "setuptools>=77", + "cython>=3.0", "scikit-image>=0.24.0", "colour-science >= 0.4.0", "build>=1.2.1", diff --git a/src/shinier/ImageProcessor.py b/src/shinier/ImageProcessor.py index 5a15e63..da83e92 100644 --- a/src/shinier/ImageProcessor.py +++ b/src/shinier/ImageProcessor.py @@ -409,7 +409,7 @@ def _compute_initial_target_histogram(self, n_bins: int = 256): target_hist = self._compute_target_hist_from_image_path(target_hist, n_bins=n_bins) else: target_hist = target_hist[:, None] if target_hist.ndim == 1 else target_hist - target_hist /= (target_hist.sum(axis=0, keepdims=True) + 1e-12) + target_hist = target_hist / (target_hist.sum(axis=0, keepdims=True) + 1e-12) if target_hist.shape[0] != n_bins: raise ValueError(f"target_hist must have {n_bins} bins, but has {target_hist.shape[0]}.") if target_hist.ndim > 1 and target_hist.shape[-1] != self.dataset.buffer.n_channels: diff --git a/src/shinier/__init__.py b/src/shinier/__init__.py index ae8c30a..2ac52e6 100644 --- a/src/shinier/__init__.py +++ b/src/shinier/__init__.py @@ -10,29 +10,53 @@ Controlling low-level image properties: The SHINE toolbox. *Behavior Research Methods, 42*(3), 671โ€“684. https://doi.org/10.3758/BRM.42.3.671 - See accompanying paper: Salvas-Hรฉbert, M., Dupuis-Roy, N., Landry, C., Charest, I. & Gosselin, F. (2025) - The SHINIER the Better: An Adaptation of the SHINE Toolbox on Python. + See accompanying paper: Salvas-Hรฉbert, M.*, Dupuis-Roy, N.*, Landry, C., Charest, I., & Gosselin, F. (2026). + SHINIER: An open-source Python package for controlling low-level image properties. + SoftwareX, 35, Article 102884. https://doi.org/10.1016/j.softx.2026.102884 """ # Metadata -__author__ = "Nicolas Dupuis-Roy" -__version__ = "0.2.0" +__author__ = "Nicolas Dupuis-Roy and Mathias Salvas-Hebert" +__version__ = "0.2.2" __email__ = "nicolas.dupuis.roy@umontreal.ca" # For direct importation from importlib import util from pathlib import Path -_HAS_CYTHON = util.find_spec("shinier._cconvolve") is not None +import sys +import warnings + +_HAS_CYTHON = False +convolve2d_direct = None +convolve2d_separable = None # This is the *package* root: src/shinier in dev, site-packages/shinier when installed DEV_ROOT = Path(__file__).resolve().parents[2] REPO_ROOT = Path(__file__).resolve().parent -if _HAS_CYTHON: - from ._cconvolve import convolve2d_direct, convolve2d_separable -else: - convolve2d_direct = None - convolve2d_separable = None +if util.find_spec("shinier._cconvolve") is not None: + try: + from ._cconvolve import convolve2d_direct, convolve2d_separable + _HAS_CYTHON = True + except Exception as exc: + try: + import numpy as _np + + numpy_version = _np.__version__ + except Exception: + numpy_version = "unavailable" + + warnings.warn( + "SHINIER could not load the optional compiled convolution extension " + "(`shinier._cconvolve`). SHINIER will keep working, but convolution-heavy " + "operations will use the slower NumPy fallback. " + f"Python: {sys.version.split()[0]}; NumPy: {numpy_version}. " + "If you want the faster compiled extension, try reinstalling SHINIER after " + "upgrading pip, setuptools, wheel, and NumPy, and make sure a C++ compiler " + f"is available. Original error: {exc!r}", + RuntimeWarning, + stacklevel=2, + ) __all__ = [ "Options", diff --git a/src/shinier/color/GamutControl.py b/src/shinier/color/GamutControl.py index 5697e28..b4df93a 100644 --- a/src/shinier/color/GamutControl.py +++ b/src/shinier/color/GamutControl.py @@ -176,7 +176,7 @@ def _reliability_from_Y01(Y01: np.ndarray, fade_width: float = 0.05) -> np.ndarr return np.clip(dist_to_edge / fade_width, 0.0, 1.0) @staticmethod - def _min_or_quantile(values: np.ndarray, quantile_threshold: float | None) -> float: + def _min_or_quantile(values: np.ndarray, quantile_threshold: Optional[float]) -> float: """Return min(values) or a lower-quantile if configured.""" if values.size == 0: return 1.0 @@ -185,7 +185,7 @@ def _min_or_quantile(values: np.ndarray, quantile_threshold: float | None) -> fl return float(np.quantile(values, quantile_threshold)) @staticmethod - def _format_quantile_msg(quantile_threshold: float | None) -> str: + def _format_quantile_msg(quantile_threshold: Optional[float]) -> str: """Format a short message describing quantile clipping.""" if quantile_threshold is None: return '' @@ -194,9 +194,9 @@ def _format_quantile_msg(quantile_threshold: float | None) -> str: def _log_image_overflow( self, kind: str, - idx: int | None, + idx: Optional[int], local_min: float, - quantile_threshold: float | None, + quantile_threshold: Optional[float], verbose: bool, ) -> None: """Log a standardized per-image overflow message. @@ -277,7 +277,7 @@ def apply_low_Y_desaturation( self, Y: np.ndarray, other: np.ndarray, - idx: int | None = None, + idx: Optional[int] = None, verbose: bool = False) -> Tuple[np.ndarray, np.ndarray]: """[S1] Optionally desaturate chroma for low-luminance pixels. diff --git a/src/shinier/color/chroma_eval.py b/src/shinier/color/chroma_eval.py index 42f1b7a..9205334 100644 --- a/src/shinier/color/chroma_eval.py +++ b/src/shinier/color/chroma_eval.py @@ -197,15 +197,13 @@ def mean_chroma_loss_pct_lab( """Compute the global relative mean chroma loss percentage in CIELAB. Both sRGB images are converted to CIE Lab, and chroma is computed for each - pixel as ``C* = sqrt(a*^2 + b*^2)``.. + pixel as ``C* = sqrt(a*^2 + b*^2)``. - The relative mean chroma loss percentage is defined as:: - - ``100 * (E[C*_before] - E[C*_after]) / E[C*_before]`` - - where ``E`` denotes the mean across all image pixels. A positive percentage - indicates a reduction in mean chroma, whereas a negative percentage indicates - an increase. + The relative mean chroma loss percentage is defined as + ``100 * (E[C*_before] - E[C*_after]) / E[C*_before]``, where ``E`` denotes + the mean across all image pixels. A positive percentage indicates a + reduction in mean chroma, whereas a negative percentage indicates an + increase. Parameters ---------- diff --git a/src/shinier/utils.py b/src/shinier/utils.py index da9d94d..a3b707b 100644 --- a/src/shinier/utils.py +++ b/src/shinier/utils.py @@ -320,16 +320,16 @@ class StimulusMasker: >>> final_mask = masker.interactive_mask(image) """ - image_size: int | tuple[int, int] + image_size: Union[int, Tuple[int, int]] cutoff_a: float - cutoff_b: float | None = None + cutoff_b: Optional[float] = None offset_a: float = 0.0 offset_b: float = 0.0 mask_type: MaskType = "feathered_disk" sigma: float = 2.0 edge_width: float = 2.0 background: float = 0.5 - output_dtype: np.dtype | type = np.float64 + output_dtype: Union[np.dtype, type] = np.float64 def mask(self) -> np.ndarray: """Generate mask as float64 in [0, 1].""" @@ -1023,7 +1023,7 @@ def get_radius_grid(x_size: int, y_size: int, legacy_mode: bool = False) -> np.n # --- polar radius, MATLAB rounding rule --- r = np.hypot(XX, YY) r_adjustment = -1 if (x_size % 2 == 1) or (y_size % 2 == 1) else 0 - r = MatlabOperators.round(r) if legacy_mode else np.round(r, decimals=0) + r_adjustment + r = (MatlabOperators.round(r) if legacy_mode else np.round(r, decimals=0)) + r_adjustment # Non-negative integer bin indices return np.clip(r, 0, None).astype(np.int64) @@ -2846,7 +2846,7 @@ def print_log(logs: List[str], log_path: Union[Path, str], log_name: Optional[st filename = Path(log_path) / log_name # Write each log to a new line in the file - with open(filename, 'w') as file: + with open(filename, 'w', encoding='utf-8') as file: for log in logs: file.write(strip_ansi(log) + '\n') @@ -5013,6 +5013,12 @@ def compute_bp2bpsim(reference: np.ndarray, enhanced: np.ndarray, n_bits: int = Values range from 0 to 1. Higher values indicate more matching bits across corresponding pixels and channels. + + References + ---------- + Rahman, H., & Paul, G. C. (2023). Tripartite sub-image histogram equalization for slightly + low contrast gray-tone image enhancement. Pattern Recognition, 134, Article 109043. + https://doi.org/10.1016/j.patcog.2022.109043 """ _check_same_shape(reference, enhanced, "compute_bp2bpsim") if n_bits < 1 or n_bits > 8: diff --git a/tests/README.md b/tests/README.md index 48ff02b..a9bb2cd 100644 --- a/tests/README.md +++ b/tests/README.md @@ -5,7 +5,7 @@ Tests are organized into **unit** and **validation (integration)** levels. --- -## โš™๏ธ Pytest Configuration +## Pytest Configuration Make sure you did to install all dev dependencies: ```pip install '.[dev]'``` ```ini @@ -26,11 +26,11 @@ All test files must either: - start with `test_`, or - end with `_test.py`\ - (e.g., `ImageDataset_test.py` โœ…) + (e.g., `ImageDataset_test.py`) --- -## ๐Ÿ“™ Markers +## Markers Use markers to select subsets of tests: @@ -38,11 +38,11 @@ Use markers to select subsets of tests: |--------------------|-----------------------------------------|-------------------------------| | `unit_tests` | Fast functional unit tests | `pytest -m unit_tests` | | `validation_tests` | Exhaustive validation (slow) | `pytest -m validation_tests` | -| `test_all_options` | Exhaustive unit tests on Options (โš  very slow; can take a few hours) | `pytest -m test_all_options` | +| `test_all_options` | Exhaustive unit tests on Options (very slow; can take a few hours) | `pytest -m test_all_options` | --- -## ๐Ÿงต Multi-Core Execution +## Multi-Core Execution Run tests in parallel automatically: @@ -58,7 +58,7 @@ pytest -n 4 -s -m unit_tests --- -## โšก Validation Tests โ€” Coverage Modes +## Validation Tests โ€” Coverage Modes `ImageProcessor_validation_test.py` supports three coverage modes selected via `COVERAGE_MODE`: @@ -88,7 +88,7 @@ Hash keys are namespaced by `COVERAGE_MODE` (e.g. `exhaustive:5` vs `pruned:5`) --- -## ๐Ÿ› ๏ธ Running Shards Locally (GNU parallel) +## Running Shards Locally (GNU parallel) Use **GNU parallel** to distribute shards across CPU cores: @@ -98,12 +98,12 @@ parallel --ungroup --jobs 8 \ pytest -s -m validation_tests' ::: 0 1 2 3 4 5 6 7 ``` -> ๐Ÿ”น `--ungroup` allows live tqdm updates in real time.\ +> `--ungroup` allows live tqdm updates in real time.\ > Without it, each shard's output is buffered until completion. --- -## ๐Ÿ–ฅ๏ธ Running on a Compute Cluster (SLURM/sbatch) +## Running on a Compute Cluster (SLURM/sbatch) Use a SLURM job array so `$SLURM_ARRAY_TASK_ID` maps directly to `SHARD_INDEX`. Create a file (e.g. `run_validation.sh`) with the following template and adapt `N_SHARDS`, @@ -146,7 +146,7 @@ START_AT=1180676 sbatch run_validation.sh --- -## ๐Ÿ”Ž Debugging Tests +## Debugging Tests ### Drop into debugger on failure @@ -168,7 +168,7 @@ pytest -m validation_tests -vv -s --tb=long --- -## ๐ŸŸก Hard vs Soft Failures in Validation Tests +## Hard vs Soft Failures in Validation Tests `ImageProcessor_validation_test.py` distinguishes two levels of failure: @@ -190,7 +190,7 @@ Soft failures produce a `.pkl` / `.json` dump in `tests/assets/tmp/` for inspect --- -## ๐Ÿ”„ Resume From a Given Combo +## Resume From a Given Combo If a bug occurs at combo 21,600 (from tqdm output): @@ -265,7 +265,7 @@ were AMBE 0.007961%, MSSIM 0.019079%, PSNR 0.014468%, BP2BPSIM 0.024168%, CI --- -## ๐Ÿค Replay a Dumped Failure +## Replay a Dumped Failure To reproduce a failed validation test: @@ -277,7 +277,7 @@ This will rebuild the same `Options`, reload selected images, and re-run the fai --- -## ๐Ÿ”น Tips +## Tips - Use `--pdb` or `--trace` for interactive debugging. - Always set `PYTHONUNBUFFERED=1` in `parallel` to force live output. @@ -293,7 +293,7 @@ This will rebuild the same `Options`, reload selected images, and re-run the fai --- -## ๐Ÿ”ง Example Workflow +## Example Workflow 1. Run sampled validation tests locally across 8 cores: @@ -314,3 +314,47 @@ This will rebuild the same `Options`, reload selected images, and re-run the fai ```bash python -m tests.tools.replay_failure path/to/failure_xxxxx.pkl ``` + +--- + +## MATLAB SHINE Comparison + +SHINIER ships with a standalone comparison tool that benchmarks the Python +implementations against the original +[MATLAB SHINE toolbox](http://www.mapageweb.umontreal.ca/gosselif/SHINE/) +across processing modes 1โ€“8. Three implementations are compared: + +| Implementation | Description | +|-----------------------|--------------------------------------------------------------------| +| `matlab_shine` | Original MATLAB SHINE toolbox, driven by a generated MATLAB script | +| `shinier_legacy` | SHINIER with `legacy_mode=True` (MATLAB-compatible behavior) | +| `shinier_modern_gray` | SHINIER defaults on grayscale (xyY luminance processing) | + +The tool produces two comparison stages: + +1. **Output comparison** โ€” pixel differences (RMSE, MAE, max abs, equal + fraction, histogram L1) between saved MATLAB and Python images. +2. **Fixed-target comparison** โ€” every implementation receives the same fixed + initial Python targets (histogram and spectrum) and each output is measured + against the target in its own processing domain. + +Requirements: a local MATLAB installation and the SHINE toolbox. + +```bash +# Complete run (asks for MATLAB/SHINE paths if not found) +bash tests/tools/run_matlab_shine_comparison.sh + +# Common overrides +MATLAB_BIN=/Applications/MATLAB_R2025a.app/bin/matlab \ +SHINE_DIR=~/toolboxes/shinetoolbox \ +MODES="2 3 4" LIMIT=8 ITERATIONS=5 \ +bash tests/tools/run_matlab_shine_comparison.sh + +# Keep all intermediate images, MATLAB scripts and .mat files +FULL_TRACKING=1 bash tests/tools/run_matlab_shine_comparison.sh +``` + +Results are written as CSV files (summary and per-image detail) under +`tmp/matlab_shine_comparison/`, and summary tables are printed to the +terminal. See the module docstring of +`tests/tools/matlab_shine_comparison.py` for the full metric definitions. diff --git a/tests/other_tests/test_install_all_pythons_on_linux.sh b/tests/other_tests/test_install_all_pythons_on_linux.sh index 5380ad3..a2757d4 100644 --- a/tests/other_tests/test_install_all_pythons_on_linux.sh +++ b/tests/other_tests/test_install_all_pythons_on_linux.sh @@ -24,6 +24,8 @@ IMAGES=( "python:3.10-slim-bookworm" "python:3.11-slim-bookworm" "python:3.12-slim-bookworm" + "python:3.13-slim-bookworm" + "python:3.14-slim-bookworm" ) ok=() diff --git a/tests/tools/matlab_shine_comparison.py b/tests/tools/matlab_shine_comparison.py new file mode 100644 index 0000000..f8184cd --- /dev/null +++ b/tests/tools/matlab_shine_comparison.py @@ -0,0 +1,1293 @@ +"""Compare SHINIER (Python) against the original MATLAB SHINE toolbox. + +Run the complete comparison with: + + bash tests/tools/run_matlab_shine_comparison.sh + +Requirements: MATLAB and the SHINE toolbox +(http://www.mapageweb.umontreal.ca/gosselif/SHINE/). + +Three implementations are compared across processing modes 1-8: + +- ``matlab_shine``: the original MATLAB SHINE toolbox (lumMatch, histMatch, + sfMatch, specMatch), driven by a generated MATLAB script that replicates + the ``SHINE.m`` composite-mode loop; +- ``shinier_legacy``: SHINIER with ``legacy_mode=True`` (MATLAB-compatible + grayscale conversion, rounding, and noise tie-breaking); +- ``shinier_modern_gray``: SHINIER defaults on grayscale (xyY luminance + processing, hybrid tie-breaking, sRGB-encoded export). + +Two comparison stages are produced: + +1. Output comparison: pixel differences between saved MATLAB and Python + images (RMSE, MAE, max abs, equal fraction, histogram L1). +2. Fixed-target comparison: every implementation receives the same fixed + initial Python targets (histogram and spectrum), and each output is + measured against the target in its own processing domain. MATLAB and + ``shinier_legacy`` share the legacy target; ``shinier_modern_gray`` uses + its own target. + +Default output is CSV-only; use ``FULL_TRACKING=1`` or ``--full-tracking`` to +keep generated inputs, PNGs, MATLAB scripts, and ``.mat`` files. The grayscale +CSV compares MATLAB's input gray image to SHINIER's internal +``ImageProcessor._initial_buffer``. PNG target metrics map saved images back +through each implementation's grayscale domain; ``png_hist_l1`` also uses +``rounded_target_hist`` and the export/import bin transform. SHINIER rows +include ``soft_clip``; MATLAB rows use SHINE's rescale/uint8 behavior. +Terminal tables use scientific notation with 3 decimals; CSVs keep full +precision. +""" + +from __future__ import annotations + +import argparse +import csv +import os +import shutil +import subprocess +import sys +import tempfile +import textwrap +from datetime import datetime +from pathlib import Path + +import numpy as np +from PIL import Image +from scipy.io import savemat + +REPO_ROOT = Path(__file__).resolve().parents[2] +DEFAULT_SHINE_DIR = Path(os.environ["SHINE_DIR"]) if os.environ.get("SHINE_DIR") else None +DEFAULT_MATLAB = Path(os.environ.get("MATLAB_BIN", "matlab")) +DEFAULT_OUTPUT_DIR = REPO_ROOT / "tmp/matlab_shine_comparison" + +if "MPLCONFIGDIR" not in os.environ: + mpl_cache = Path(tempfile.gettempdir()) / "shinier_matplotlib" + mpl_cache.mkdir(parents=True, exist_ok=True) + os.environ["MPLCONFIGDIR"] = str(mpl_cache) + +sys.path.insert(0, str(REPO_ROOT / "src")) + +from shinier import ImageDataset, ImageProcessor, Options # noqa: E402 +from shinier.color.Converter import ColorConverter, rgb2gray # noqa: E402 +from shinier.utils import ( # noqa: E402 + _crop_after_fft, + compute_rmse, + compute_tvd_hist, + get_radius_grid, + image_spectrum, + MatlabOperators, + pol2cart, + rounded_target_hist, + rotational_avg, + soft_clip, + uint8_plus, +) + +TARGET_REFERENCE = ( + "Initial Python target. MATLAB and SHINIER legacy use the legacy target; " + "SHINIER modern-gray uses its own modern-gray target. Composite modes are " + "evaluated against this fixed initial target. In this target-injection run, " + "MATLAB and Python are given the fixed target; normal composite processing " + "without injected targets may use moving targets internally." +) +ANSI_BOLD_RED = "\033[1;31m" +ANSI_RESET = "\033[0m" +MATLAB_IMPLEMENTATION = "matlab_shine" +LEGACY_IMPLEMENTATION = "shinier_legacy" +MODERN_GRAY_IMPLEMENTATION = "shinier_modern_gray" +PYTHON_IMPLEMENTATIONS = { + LEGACY_IMPLEMENTATION: ("legacy", True), + MODERN_GRAY_IMPLEMENTATION: ("modern_gray", False), +} + +HIST_MODES = {2, 5, 6, 7, 8} +SF_MODES = {3, 5, 7} +SPECTRUM_MODES = {4, 6, 8} +TARGET_MODES = HIST_MODES | SF_MODES | SPECTRUM_MODES + +TARGET_METRIC_LABELS = { + "internal_buffer_hist_l1_to_python_target": "int_hist_l1", + "exported_png_hist_l1_to_python_target": "png_hist_l1", + "internal_buffer_sf_rmse_to_python_target": "int_sf_rmse", + "exported_png_sf_rmse_to_python_target": "png_sf_rmse", + "pre_range_fourier_spectrum_rmse_to_python_target": "fft_spec", + "pre_range_image_spectrum_rmse_to_python_target": "pre_spec", + "pre_range_out_of_range_fraction": "pre_oor", + "post_soft_clip_spectrum_rmse_to_python_target": "clip_spec", + "internal_buffer_spectrum_rmse_to_python_target": "int_spec_rmse", + "exported_png_spectrum_rmse_to_python_target": "png_spec_rmse", +} +PIXEL_COLUMNS = "images mean_rmse mean_mae max_abs mean_equal_fraction mean_hist_l1".split() +OUTPUT_COMPARISON_COLUMNS = "mode reference python_output".split() + PIXEL_COLUMNS[1:] +TARGET_TABLE_COLUMNS = ["mode", "implementation", "target"] + [ + f"mean_{key}" for key in TARGET_METRIC_LABELS +] +TARGET_TABLE_LABELS = {f"mean_{key}": label for key, label in TARGET_METRIC_LABELS.items()} +# CSV-only keys extend the displayed metrics with the signed out-of-range fractions. +TARGET_METRIC_KEYS = ( + *TARGET_METRIC_LABELS, + "pre_range_below_zero_fraction", + "pre_range_above_one_fraction", +) +TARGET_TABLE_NOTES = [ + "Legend:", + " int_* = metric on SHINIER ImageProcessor._final_buffer; blank for MATLAB.", + " png_* = exported PNG mapped back to the implementation target domain.", + " fft/pre/clip_spec = SHINIER-only spectrum checks before/after soft_clip; pre_oor = pre-correction out-of-range fraction.", + " clip_spec vs int_spec_rmse gap = final rescaling for modes 4/6, final hist_match for mode 8.", + " SHINIER uses soft_clip, which lowers hard clipping but can worsen fixed-target RMSE; MATLAB uses SHINE rescale/uint8.", + " png_hist_l1 uses rounded_target_hist plus the same export/import bin transform.", +] +TARGET_WARNING_LINES = [ + "target metrics use fixed initial Python targets; normal composite runs may use moving targets internally.", + "MATLAB/legacy share the legacy target; modern-gray uses its own target.", + "SHINIER rows include soft_clip; MATLAB rows use SHINE rescale/uint8.", + "png_hist_l1 uses rounded_target_hist plus the export/import bin mapping.", +] +RUN_INFO_KEYS = [ + "run_root", + "matlab_runner", + "modes", + "images", + "use_python_targets", + "skip_matlab", + "prepare_only", +] +INPUT_GRAY_TITLE = "MATLAB input grayscale vs SHINIER internal initial buffer" +TargetSet = dict[str, np.ndarray | int | str | bool] +TargetSets = dict[str, TargetSet] + + +def parse_args() -> argparse.Namespace: + parser = argparse.ArgumentParser( + description=__doc__, + formatter_class=argparse.RawDescriptionHelpFormatter, + ) + add = parser.add_argument + add("--input-dir", type=Path, default=REPO_ROOT / "tests/assets/SAMPLE_64X64") + add("--output-dir", type=Path, default=DEFAULT_OUTPUT_DIR) + add("--run-root", type=Path) + add("--shine-dir", type=Path, default=DEFAULT_SHINE_DIR) + add("--matlab-bin", type=Path, default=DEFAULT_MATLAB) + add("--modes", type=int, nargs="+", default=list(range(1, 9))) + add("--iterations", type=int, default=5) + add("--seed", type=int, default=42) + add("--limit", type=int, default=8) + for flag in ( + "--skip-matlab", + "--prepare-only", + "--use-python-targets", + "--quiet-run-info", + "--full-tracking", + ): + add(flag, action="store_true") + return parser.parse_args() + + +def quote_matlab_path(path: Path) -> str: + return str(path).replace("'", "''") + + +def prepare_inputs(source_dir: Path, output_root: Path, limit: int) -> list[Path]: + input_dir = output_root / "input" + input_dir.mkdir(parents=True, exist_ok=True) + paths = sorted(source_dir.glob("*.png"))[:limit] + if not paths: + raise FileNotFoundError(f"No PNG inputs found in {source_dir}") + for path in paths: + shutil.copy2(path, input_dir / path.name) + return [input_dir / path.name for path in paths] + + +def write_matlab_runner( + *, + output_root: Path, + input_dir: Path, + shine_dir: Path | None, + modes: list[int], + iterations: int, + seed: int, + use_python_targets: bool = False, +) -> Path: + runner = output_root / "run_matlab_shine.m" + modes_text = " ".join(str(mode) for mode in modes) + shine_addpath = f"addpath('{quote_matlab_path(shine_dir)}');" if shine_dir else "" + shine_check = ( + "required_shine = {'lumMatch', 'histMatch', 'sfMatch', 'specMatch'};\n" + "missing_shine = {};\n" + "for si = 1:numel(required_shine)\n" + " if exist(required_shine{si}, 'file') ~= 2\n" + " missing_shine{end + 1} = required_shine{si}; %#ok\n" + " end\n" + "end\n" + "if ~isempty(missing_shine)\n" + " error(['Missing SHINE functions: ', strjoin(missing_shine, ', '), ...\n" + " '. Set SHINE_DIR to the SHINE toolbox folder or download SHINE from ', ...\n" + " 'http://www.mapageweb.umontreal.ca/gosselif/SHINE/']);\n" + "end" + ) + target_load = "" + hist_arg = "" + spectrum_arg = "" + if use_python_targets: + target_load = ( + f"targets = load('{quote_matlab_path(output_root / 'python_targets.mat')}');" + ) + hist_arg = ", targets.target_hist_counts" + spectrum_arg = ", targets.target_spectrum" + runner.write_text( + textwrap.dedent( + f""" + {shine_addpath} + {shine_check} + input_dir = '{quote_matlab_path(input_dir)}'; + output_root = '{quote_matlab_path(output_root / "matlab")}'; + modes = [{modes_text}]; + {target_load} + files = dir(fullfile(input_dir, '*.png')); + [~, order] = sort({{files.name}}); + files = files(order); + + original = cell(1, numel(files)); + input_gray_dir = fullfile(output_root, 'input_gray'); + if ~exist(input_gray_dir, 'dir') + mkdir(input_gray_dir); + end + for k = 1:numel(files) + im = imread(fullfile(input_dir, files(k).name)); + if ndims(im) == 3 + im = rgb2gray(im); + end + original{{k}} = im; + imwrite(uint8(im), fullfile(input_gray_dir, files(k).name)); + end + + hist_corr_summary = {{}}; + for mi = 1:numel(modes) + mode = modes(mi); + rand('seed', {seed}); + images = original; + if mode >= 5 + n_iter = {iterations}; + else + n_iter = 1; + end + for iter = 1:n_iter + switch mode + case 1 + images = lumMatch(images); + case {{2, 5, 6}} + [shine_log, images] = evalc('histMatch(images, 0{hist_arg})'); + hist_corr_summary = append_hist_corr( ... + hist_corr_summary, mode, iter, 'hist_before_spectrum', shine_log); + end + switch mode + case {{3, 5, 7}} + images = sfMatch(images, 1{spectrum_arg}); + case {{4, 6, 8}} + images = specMatch(images, 1{spectrum_arg}); + end + switch mode + case {{7, 8}} + [shine_log, images] = evalc('histMatch(images, 0{hist_arg})'); + hist_corr_summary = append_hist_corr( ... + hist_corr_summary, mode, iter, 'hist_after_spectrum', shine_log); + end + end + + out_dir = fullfile(output_root, sprintf('mode_%d', mode)); + if ~exist(out_dir, 'dir') + mkdir(out_dir); + end + for k = 1:numel(files) + imwrite(uint8(images{{k}}), fullfile(out_dir, files(k).name)); + end + end + + print_hist_corr_summary(hist_corr_summary); + + function rows = append_hist_corr(rows, mode, iter, stage, shine_log) + tokens = regexp( ... + shine_log, ... + 'Correlation between processed image and target histogram:\\s*([0-9.eE+-]+)', ... + 'tokens'); + if isempty(tokens) + return; + end + values = zeros(1, numel(tokens)); + for ti = 1:numel(tokens) + values(ti) = str2double(tokens{{ti}}{{1}}); + end + rows(end + 1, :) = {{ ... + mode, iter, stage, numel(values), min(values), mean(values), max(values) ... + }}; + end + + function print_hist_corr_summary(rows) + if isempty(rows) + return; + end + fprintf('\\nMATLAB SHINE histMatch correlation summary\\n'); + fprintf('------+-----------+----------------------+--------+-----------+-----------+-----------\\n'); + fprintf('mode | iteration | stage | images | min_corr | mean_corr | max_corr\\n'); + fprintf('------+-----------+----------------------+--------+-----------+-----------+-----------\\n'); + for ri = 1:size(rows, 1) + fprintf( ... + '%-5d | %-9d | %-20s | %-6d | %.3e | %.3e | %.3e\\n', ... + rows{{ri, 1}}, rows{{ri, 2}}, rows{{ri, 3}}, rows{{ri, 4}}, ... + rows{{ri, 5}}, rows{{ri, 6}}, rows{{ri, 7}}); + end + fprintf('------+-----------+----------------------+--------+-----------+-----------+-----------\\n'); + end + """ + ).strip() + + "\n" + ) + return runner + + +def _build_options( + *, + input_dir: Path, + output_dir: Path, + mode: int, + iterations: int, + seed: int, + target_hist: np.ndarray | None = None, + target_spectrum: np.ndarray | None = None, + legacy_mode: bool = True, +) -> Options: + return Options( + input_folder=input_dir, + output_folder=output_dir, + mode=mode, + legacy_mode=legacy_mode, + as_gray=True, + # Keep False here: True uses equal-channel grayscale; False uses the + # rec601 color-treatment path (legacy rgb2gray or modern xyY/Y). + linear_luminance=False, + rec_standard=1, + iterations=iterations, + seed=seed, + target_hist=target_hist, + target_spectrum=target_spectrum, + verbose=-1, + ) + + +def compute_python_targets( + *, + input_dir: Path, + output_root: Path, + seed: int, + legacy_mode: bool = True, + target_name: str = "legacy", + export_matlab_targets: bool = True, +) -> TargetSet: + target_dir = output_root / f"python_target_source_{target_name}" + target_dir.mkdir(parents=True, exist_ok=True) + # Only the initial targets are read from this run; they are computed before + # the iteration loop, so a single iteration is enough. + options = _build_options( + input_dir=input_dir, + output_dir=target_dir, + mode=8, + iterations=1, + seed=seed, + legacy_mode=legacy_mode, + ) + proc = ImageProcessor(dataset=ImageDataset(options=options), options=options) + + target_hist = np.asarray(proc._initial_targets["hist"], dtype=np.float64) + target_hist_1d = target_hist[:, 0] if target_hist.ndim == 2 else target_hist + first_image = read_saved_gray(sorted(input_dir.glob("*.png"))[0]) + n_pixels = int(first_image.size) + target_hist_freq = rounded_target_hist(target_hist_1d, n_pixels) + target_hist_counts = np.rint(target_hist_freq * n_pixels).astype(np.int64) + + target_spectrum = np.asarray(proc._initial_targets["spectrum"], dtype=np.float64).squeeze() + radius = get_radius_grid( + target_spectrum.shape[0], + target_spectrum.shape[1], + legacy_mode=legacy_mode, + ) + target_sf = rotational_avg(target_spectrum, radius) + + if export_matlab_targets: + target_path = output_root / "python_targets.mat" + savemat( + target_path, + { + "target_hist_freq": target_hist_freq[:, None], + "target_hist_counts": target_hist_counts[:, None], + "target_spectrum": target_spectrum, + "target_sf": target_sf[:, None], + "n_pixels": np.array([[n_pixels]], dtype=np.int64), + }, + ) + (output_root / "python_target_reference.txt").write_text(TARGET_REFERENCE + "\n") + return { + "target_name": target_name, + "legacy_mode": legacy_mode, + "target_hist_freq": target_hist_freq, + "target_hist_counts": target_hist_counts, + "target_spectrum": target_spectrum, + "target_sf": target_sf, + "n_pixels": n_pixels, + } + + +def run_matlab(matlab_bin: Path, runner: Path) -> None: + cmd = [str(matlab_bin), "-batch", f"run('{quote_matlab_path(runner)}')"] + completed = subprocess.run(cmd, text=True, capture_output=True) + if completed.stdout: + print(completed.stdout) + if completed.stderr: + print(completed.stderr, file=sys.stderr) + if completed.returncode != 0: + print( + f"MATLAB exited with status {completed.returncode}; " + "continuing so output checks can decide.", + file=sys.stderr, + ) + + +def run_python_processor( + *, + input_dir: Path, + output_root: Path, + implementation: str, + mode: int, + iterations: int, + seed: int, + target_hist: np.ndarray | None = None, + target_spectrum: np.ndarray | None = None, + keep_internal: bool = False, +) -> tuple[Path, ImageProcessor]: + out_dir = output_root / "python" / implementation / f"mode_{mode}" + out_dir.mkdir(parents=True, exist_ok=True) + options = _build_options( + input_dir=input_dir, + output_dir=out_dir, + mode=mode, + iterations=iterations, + seed=seed, + target_hist=target_hist, + target_spectrum=target_spectrum, + legacy_mode=PYTHON_IMPLEMENTATIONS[implementation][1], + ) + processor = ImageProcessor( + dataset=ImageDataset(options=options), + options=options, + from_unit_test=keep_internal, + ) + if keep_internal: + processor.process() + processor.print_log_results() + if not getattr(processor.dataset.images, "has_list_array", False): + processor.dataset.save_images() + return out_dir, processor + + +def read_saved_gray(path: Path, implementation: str = MATLAB_IMPLEMENTATION) -> np.ndarray: + image = np.asarray(Image.open(path)) + if image.ndim == 2: + return image.astype(np.float64) + + rgb = image[..., :3] + if implementation in {MATLAB_IMPLEMENTATION, LEGACY_IMPLEMENTATION}: + return MatlabOperators.uint8( + rgb2gray(rgb, weighting_standard="rec601", matlab_601=True) + ).astype(np.float64) + + if implementation == MODERN_GRAY_IMPLEMENTATION: + if np.array_equal(rgb[..., 0], rgb[..., 1]) and np.array_equal(rgb[..., 0], rgb[..., 2]): + return rgb[..., 0].astype(np.float64) + return np.rint( + np.clip(rgb2gray(rgb, weighting_standard="rec601", matlab_601=False), 0, 255) + ).astype(np.float64) + + raise ValueError(f"Unknown saved-output grayscale conversion for: {implementation}") + + +def read_saved_target_image(path: Path, implementation: str = MATLAB_IMPLEMENTATION) -> np.ndarray: + if implementation != MODERN_GRAY_IMPLEMENTATION: + return read_saved_gray(path, implementation=implementation) + + image = np.asarray(Image.open(path)) + if image.ndim == 2: + image = np.dstack([image, image, image]) + rgb = image[..., :3].astype(np.float64) / 255.0 + converter = ColorConverter(rec_standard="rec601") + return converter.sRGB_to_xyY(rgb)[..., 2] * 255.0 + + +def load_shinier_initial_buffer( + *, + input_dir: Path, + output_root: Path, + seed: int, + full_tracking: bool, +) -> dict[str, np.ndarray]: + probe_dir = output_root / "python" / "internal_input_probe" + probe_dir.mkdir(parents=True, exist_ok=True) + options = _build_options( + input_dir=input_dir, + output_dir=probe_dir, + mode=1, + iterations=1, + seed=seed, + ) + processor = ImageProcessor(dataset=ImageDataset(options=options), options=options) + input_paths = sorted(input_dir.glob("*.png")) + internal = { + path.name: np.asarray(processor._initial_buffer[idx], dtype=np.float64) + for idx, path in enumerate(input_paths) + } + processor.dataset.close() + + if full_tracking: + out_dir = output_root / "python" / "internal_input_gray" + out_dir.mkdir(parents=True, exist_ok=True) + for name, image in internal.items(): + Image.fromarray(np.clip(np.rint(image), 0, 255).astype(np.uint8)).save(out_dir / name) + return internal + + +def hist_l1(a: np.ndarray, b: np.ndarray) -> float: + """L1 distance between two image histograms (= 2 x total variation distance).""" + ha = np.bincount(a.astype(np.uint8).ravel(), minlength=256).astype(np.float64) + hb = np.bincount(b.astype(np.uint8).ravel(), minlength=256).astype(np.float64) + ha /= max(ha.sum(), 1.0) + hb /= max(hb.sum(), 1.0) + return 2.0 * compute_tvd_hist(ha, hb) + + +def hist_l1_to_target(image: np.ndarray, target_hist_freq: np.ndarray) -> float: + """L1 distance between an image histogram and a target distribution.""" + binned = np.clip(np.rint(image), 0, 255).astype(np.uint8) + hist = np.bincount(binned.ravel(), minlength=256).astype(np.float64) + hist /= max(hist.sum(), 1.0) + target = np.asarray(target_hist_freq, dtype=np.float64).reshape(-1) + target = target / target.sum() + return 2.0 * compute_tvd_hist(hist, target) + + +def exported_hist_target( + target_hist_freq: np.ndarray, + *, + implementation: str, + n_pixels: int | None = None, +) -> np.ndarray: + target = np.asarray(target_hist_freq, dtype=np.float64).reshape(-1) + target = target / target.sum() + if n_pixels is not None: + target = rounded_target_hist(target, int(n_pixels)) + + if implementation != MODERN_GRAY_IMPLEMENTATION: + return target + + levels = np.arange(256, dtype=np.float64) + converter = ColorConverter(rec_standard="rec601") + srgb = converter.linRGB_to_sRGB( + np.dstack([levels / 255, levels / 255, levels / 255]) + )[..., 0] * 255 + saved = uint8_plus(srgb) + target_domain = converter.sRGB_to_xyY( + np.dstack([saved, saved, saved]).astype(np.float64) / 255 + )[..., 2] * 255 + mapped_bins = np.clip(np.rint(target_domain), 0, 255).astype(np.uint8).ravel() + + mapped = np.zeros_like(target) + for src_bin, dst_bin in enumerate(mapped_bins): + mapped[int(dst_bin)] += target[src_bin] + mapped /= mapped.sum() + return mapped + + +def spectrum_rmse_to_target(image: np.ndarray, target_spectrum: np.ndarray) -> float: + """RMSE between an image's Fourier magnitude spectrum and a target spectrum.""" + magnitude, _ = image_spectrum(image / 255.0, rescale=False) + observed = magnitude.squeeze().astype(np.float64) + target = np.asarray(target_spectrum, dtype=np.float64).squeeze() + return float(compute_rmse(observed, target)) + + +def pre_range_spectrum_metrics( + *, + processor: ImageProcessor, + input_dir: Path, + mode: int, + targets: TargetSet, +) -> dict[str, dict[str, float | str]]: + if mode not in SPECTRUM_MODES: + return {} + + target = np.asarray(targets["target_spectrum"], dtype=np.float64) + if target.ndim == 2: + target = target[..., None] + + rows: dict[str, dict[str, float | str]] = {} + for idx, input_path in enumerate(sorted(input_dir.glob("*.png"))): + phase = np.asarray(processor.dataset.phases[idx], dtype=np.float64) + if phase.ndim == 2: + phase = phase[..., None] + + original_shape = np.asarray(processor._initial_buffer[idx]).shape[:2] + pre_range_channels: list[np.ndarray] = [] + fourier_rmses: list[float] = [] + for channel in range(target.shape[-1]): + xx, yy = pol2cart(target[:, :, channel], phase[:, :, channel]) + reconstructed = np.real(np.fft.ifft2(np.fft.ifftshift(xx + yy * 1j))) + reconstructed_mag = np.abs(np.fft.fftshift(np.fft.fft2(reconstructed))) + fourier_rmses.append(float(compute_rmse(reconstructed_mag, target[:, :, channel]))) + pre_range_channels.append(_crop_after_fft(reconstructed, original_shape)) + + pre_range01 = np.stack(pre_range_channels, axis=-1).squeeze() + pre_range255 = pre_range01 * 255.0 + try: + pre_range_spec_rmse: float | str = spectrum_rmse_to_target(pre_range255, target) + except ValueError: + pre_range_spec_rmse = "" + + if pre_range01.min() < 0.0 or pre_range01.max() > 1.0: + post_clip01 = soft_clip( + pre_range01, + min_value=0.0, + max_value=1.0, + max_percent=0.01, + verbose=False, + ) + else: + post_clip01 = pre_range01 + post_clip_spec_rmse = spectrum_rmse_to_target(post_clip01 * 255.0, target) + + rows[input_path.name] = { + "pre_range_fourier_spectrum_rmse_to_python_target": float(np.mean(fourier_rmses)), + "pre_range_image_spectrum_rmse_to_python_target": pre_range_spec_rmse, + "post_soft_clip_spectrum_rmse_to_python_target": post_clip_spec_rmse, + "pre_range_below_zero_fraction": float(np.mean(pre_range01 < 0.0)), + "pre_range_above_one_fraction": float(np.mean(pre_range01 > 1.0)), + "pre_range_out_of_range_fraction": float( + np.mean((pre_range01 < 0.0) | (pre_range01 > 1.0)) + ), + } + return rows + + +def sf_rmse_to_target(image: np.ndarray, target_sf: np.ndarray, legacy_mode: bool = True) -> float: + """RMSE between an image's rotational-average SF profile and a target profile.""" + magnitude, _ = image_spectrum(image / 255.0, rescale=False) + observed_mag = magnitude.squeeze().astype(np.float64) + radius = get_radius_grid(observed_mag.shape[0], observed_mag.shape[1], legacy_mode=legacy_mode) + observed_sf = rotational_avg(observed_mag, radius) + target = np.asarray(target_sf, dtype=np.float64).reshape(-1) + n = min(observed_sf.size, target.size) + return float(compute_rmse(observed_sf[:n], target[:n])) + + +def pixel_difference_metrics(a: np.ndarray, b: np.ndarray) -> dict[str, object]: + diff = a - b + abs_diff = np.abs(diff) + return { + "rmse": float(compute_rmse(a, b)), + "mae": float(np.mean(abs_diff)), + "max_abs": float(abs_diff.max()), + "nonzero_pixels": int(np.count_nonzero(abs_diff)), + "equal_fraction": float(np.mean(abs_diff == 0)), + "hist_l1": hist_l1(a, b), + } + + +def mean_value(rows: list[dict[str, object]], key: str) -> float: + return float(np.mean([float(row[key]) for row in rows])) + + +def group_rows(rows: list[dict[str, object]], *keys: str): + groups: dict[tuple[object, ...], list[dict[str, object]]] = {} + for row in rows: + groups.setdefault(tuple(row.get(key, "") for key in keys), []).append(row) + return sorted(groups.items()) + + +def summarize_pixel_rows(rows: list[dict[str, object]]) -> dict[str, object]: + return { + "images": len(rows), + "mean_rmse": mean_value(rows, "rmse"), + "mean_mae": mean_value(rows, "mae"), + "max_abs": float(max(float(row["max_abs"]) for row in rows)), + "mean_equal_fraction": mean_value(rows, "equal_fraction"), + "mean_hist_l1": mean_value(rows, "hist_l1"), + } + + +def target_row( + *, + mode: int, + implementation: str, + targets: TargetSet, + image_name: str, + metrics: dict[str, object], +) -> dict[str, object]: + return { + "mode": mode, + "implementation": implementation, + "target": str(targets.get("target_name", "legacy")), + "image": image_name, + **{key: "" for key in TARGET_METRIC_KEYS}, + **metrics, + } + + +def target_metrics( + image: np.ndarray, + *, + mode: int, + targets: TargetSet, + prefix: str, + hist_target_freq: np.ndarray | None = None, +) -> dict[str, object]: + metrics: dict[str, object] = {} + if mode in HIST_MODES: + metrics[f"{prefix}_hist_l1_to_python_target"] = hist_l1_to_target( + image, + hist_target_freq if hist_target_freq is not None else targets["target_hist_freq"], + ) + if mode in SF_MODES: + metrics[f"{prefix}_sf_rmse_to_python_target"] = sf_rmse_to_target( + image, + targets["target_sf"], + legacy_mode=bool(targets.get("legacy_mode", True)), + ) + if mode in SPECTRUM_MODES: + metrics[f"{prefix}_spectrum_rmse_to_python_target"] = spectrum_rmse_to_target( + image, targets["target_spectrum"] + ) + return metrics + + +def exported_target_metrics( + image: np.ndarray, + *, + implementation: str, + mode: int, + targets: TargetSet, +) -> dict[str, object]: + hist_target_freq = None + if mode in HIST_MODES: + hist_target_freq = exported_hist_target( + targets["target_hist_freq"], + implementation=implementation, + n_pixels=int(targets["n_pixels"]), + ) + return target_metrics( + image, + mode=mode, + targets=targets, + prefix="exported_png", + hist_target_freq=hist_target_freq, + ) + + +def compare_dirs( + matlab_dir: Path, + python_dir: Path, + mode: int, + implementation: str, +) -> list[dict[str, object]]: + rows: list[dict[str, object]] = [] + for matlab_path in sorted(matlab_dir.glob("*.png")): + python_path = python_dir / matlab_path.name + if not python_path.exists(): + raise FileNotFoundError(f"Missing Python output: {python_path}") + m = read_saved_gray(matlab_path, implementation=MATLAB_IMPLEMENTATION) + p = read_saved_gray(python_path, implementation=implementation) + rows.append({ + "mode": mode, + "reference": MATLAB_IMPLEMENTATION, + "python_output": implementation, + "image": matlab_path.name, + **pixel_difference_metrics(p, m), + }) + return rows + + +def compare_input_grayscale( + input_dir: Path, + output_root: Path, + *, + seed: int, + full_tracking: bool = False, +) -> list[dict[str, object]]: + matlab_gray_dir = output_root / "matlab" / "input_gray" + if not matlab_gray_dir.exists(): + return [] + shinier_internal = load_shinier_initial_buffer( + input_dir=input_dir, + output_root=output_root, + seed=seed, + full_tracking=full_tracking, + ) + + rows: list[dict[str, object]] = [] + for matlab_path in sorted(matlab_gray_dir.glob("*.png")): + if matlab_path.name not in shinier_internal: + raise FileNotFoundError( + f"Missing SHINIER internal input buffer for {matlab_path.name}" + ) + m = read_saved_gray(matlab_path, implementation=MATLAB_IMPLEMENTATION) + p = shinier_internal[matlab_path.name] + rows.append({"image": matlab_path.name, **pixel_difference_metrics(p, m)}) + write_csv(output_root / "input_grayscale_comparison.csv", rows) + return rows + + +def aggregate(rows: list[dict[str, object]]) -> list[dict[str, object]]: + out = [] + for (mode, reference, python_output), items in group_rows( + rows, "mode", "reference", "python_output" + ): + summary = summarize_pixel_rows(items) + out.append( + { + "mode": int(mode), + "reference": reference, + "python_output": python_output, + **{key: summary[key] for key in PIXEL_COLUMNS if key != "images"}, + } + ) + return out + + +def compare_outputs_to_python_targets( + *, + output_dir: Path, + implementation: str, + mode: int, + targets: TargetSet, +) -> list[dict[str, object]]: + if mode not in TARGET_MODES: + return [] + rows: list[dict[str, object]] = [] + for image_path in sorted(output_dir.glob("*.png")): + image = read_saved_target_image(image_path, implementation=implementation) + rows.append( + target_row( + mode=mode, + implementation=implementation, + targets=targets, + image_name=image_path.name, + metrics=exported_target_metrics( + image, + implementation=implementation, + mode=mode, + targets=targets, + ), + ) + ) + return rows + + +def compare_processor_to_python_targets( + *, + processor: ImageProcessor, + input_dir: Path, + output_dir: Path, + implementation: str, + mode: int, + targets: TargetSet, +) -> list[dict[str, object]]: + if mode not in TARGET_MODES: + return [] + rows: list[dict[str, object]] = [] + pre_range_by_image = pre_range_spectrum_metrics( + processor=processor, + input_dir=input_dir, + mode=mode, + targets=targets, + ) + for idx, input_path in enumerate(sorted(input_dir.glob("*.png"))): + image = np.asarray(processor._final_buffer[idx]).squeeze().astype(np.float64) + metrics = { + **target_metrics(image, mode=mode, targets=targets, prefix="internal_buffer"), + **pre_range_by_image.get(input_path.name, {}), + } + + saved_path = output_dir / input_path.name + if saved_path.exists(): + exported_image = read_saved_target_image(saved_path, implementation=implementation) + metrics.update( + exported_target_metrics( + exported_image, + implementation=implementation, + mode=mode, + targets=targets, + ) + ) + rows.append( + target_row( + mode=mode, + implementation=implementation, + targets=targets, + image_name=input_path.name, + metrics=metrics, + ) + ) + return rows + + +def aggregate_target_rows(rows: list[dict[str, object]]) -> list[dict[str, object]]: + out: list[dict[str, object]] = [] + for (mode, implementation, target), items in group_rows( + rows, "mode", "implementation", "target" + ): + summary: dict[str, object] = { + "mode": int(mode), + "implementation": implementation, + "target": target, + } + for key in TARGET_METRIC_KEYS: + vals = [float(item[key]) for item in items if item[key] != ""] + summary[f"mean_{key}"] = float(np.mean(vals)) if vals else "" + out.append(summary) + return out + + +def write_csv(path: Path, rows: list[dict[str, object]]) -> None: + if not rows: + return + path.parent.mkdir(parents=True, exist_ok=True) + with path.open("w", newline="", encoding="utf-8") as f: + writer = csv.DictWriter(f, fieldnames=list(rows[0].keys())) + writer.writeheader() + writer.writerows(rows) + + +def cleanup_tracking_artifacts(run_root: Path) -> None: + if not run_root.exists(): + return + for path in run_root.iterdir(): + if path.is_dir(): + shutil.rmtree(path) + elif path.suffix.lower() != ".csv": + path.unlink() + + +INTEGER_COLUMNS = {"mode", "images", "nonzero_pixels"} + + +def format_table_value(column: str, value: object) -> str: + if value == "": + return "-" + if isinstance(value, bool): + return str(value).lower() + if isinstance(value, (int, np.integer)): + return str(int(value)) + if column in INTEGER_COLUMNS and isinstance(value, (int, float, np.integer, np.floating)): + return str(int(value)) + if isinstance(value, (int, float, np.integer, np.floating)): + return f"{float(value):.3e}" + return str(value) + + +def print_table( + title: str, + rows: list[dict[str, object]], + columns: list[str], + *, + column_labels: dict[str, str] | None = None, + notes: list[str] | None = None, +) -> None: + if not rows: + return + labels = column_labels or {} + headers = [labels.get(column, column) for column in columns] + formatted_rows = [ + [format_table_value(column, row.get(column, "")) for column in columns] + for row in rows + ] + widths = [ + max(len(headers[index]), *(len(row[index]) for row in formatted_rows)) + for index, column in enumerate(columns) + ] + header = " | ".join(headers[index].ljust(widths[index]) for index, _ in enumerate(columns)) + separator = "-+-".join("-" * width for width in widths) + + print() + print(title) + print(separator) + print(header) + print(separator) + for row in formatted_rows: + print(" | ".join(value.ljust(widths[index]) for index, value in enumerate(row))) + print(separator) + if notes: + for note in notes: + print(note) + + +def print_key_values(title: str, row: dict[str, object], keys: list[str]) -> None: + key_width = max(len(key) for key in keys) + + print() + print(title) + print("-" * len(title)) + for key in keys: + value = row.get(key, "") + if isinstance(value, bool): + text = str(value).lower() + else: + text = str(value) + print(f"{key.ljust(key_width)} : {text}") + + +def prepare_run(args: argparse.Namespace) -> tuple[Path, Path, Path]: + run_root = ( + args.run_root or args.output_dir / datetime.now().strftime("%Y%m%d_%H%M%S") + ).resolve() + run_root.mkdir(parents=True, exist_ok=True) + input_dir = run_root / "input" + if not input_dir.exists() or not list(input_dir.glob("*.png")): + copied_inputs = prepare_inputs(args.input_dir, run_root, args.limit) + input_dir = copied_inputs[0].parent + + runner = write_matlab_runner( + output_root=run_root, + input_dir=input_dir, + shine_dir=args.shine_dir, + modes=args.modes, + iterations=args.iterations, + seed=args.seed, + use_python_targets=args.use_python_targets, + ) + if not args.quiet_run_info: + print_key_values( + "Prepared diagnostic run", + { + "run_root": run_root, + "matlab_runner": runner, + "modes": " ".join(str(mode) for mode in args.modes), + "images": len(list(input_dir.glob("*.png"))), + "use_python_targets": args.use_python_targets, + "skip_matlab": args.skip_matlab, + "prepare_only": args.prepare_only, + }, + RUN_INFO_KEYS, + ) + return run_root, input_dir, runner + + +def compute_target_sets( + *, + input_dir: Path, + run_root: Path, + seed: int, +) -> TargetSets: + configs = [ + ("legacy", True, True), + ("modern_gray", False, False), + ] + return { + name: compute_python_targets( + input_dir=input_dir, + output_root=run_root, + seed=seed, + legacy_mode=legacy_mode, + target_name=name, + export_matlab_targets=export_matlab_targets, + ) + for name, legacy_mode, export_matlab_targets in configs + } + + +def print_input_grayscale_results(rows: list[dict[str, object]], run_root: Path) -> None: + if rows: + print(f"input_gray_summary={run_root / 'input_grayscale_comparison.csv'}") + print_table(INPUT_GRAY_TITLE, [summarize_pixel_rows(rows)], PIXEL_COLUMNS) + return + print_table( + INPUT_GRAY_TITLE, + [ + { + "status": "not_available", + "reason": "MATLAB input_gray outputs were not found", + } + ], + ["status", "reason"], + ) + + +def matlab_target_rows( + *, + run_root: Path, + mode: int, + skip_matlab: bool, + target_sets: TargetSets, +) -> list[dict[str, object]]: + matlab_dir = run_root / "matlab" / f"mode_{mode}" + if skip_matlab and not matlab_dir.exists(): + return [] + if not matlab_dir.exists(): + raise FileNotFoundError(f"Missing MATLAB outputs for mode {mode}: {matlab_dir}") + return compare_outputs_to_python_targets( + output_dir=matlab_dir, + implementation=MATLAB_IMPLEMENTATION, + mode=mode, + targets=target_sets["legacy"], + ) + + +def python_target_rows( + *, + args: argparse.Namespace, + input_dir: Path, + run_root: Path, + mode: int, + target_sets: TargetSets, +) -> list[dict[str, object]]: + rows: list[dict[str, object]] = [] + for implementation in PYTHON_IMPLEMENTATIONS: + targets = target_sets[PYTHON_IMPLEMENTATIONS[implementation][0]] + python_dir, processor = run_python_processor( + input_dir=input_dir, + output_root=run_root, + implementation=implementation, + mode=mode, + iterations=args.iterations, + seed=args.seed, + target_hist=np.asarray(targets["target_hist_freq"], dtype=np.float64), + target_spectrum=np.asarray(targets["target_spectrum"], dtype=np.float64), + keep_internal=True, + ) + try: + rows.extend( + compare_processor_to_python_targets( + processor=processor, + input_dir=input_dir, + output_dir=python_dir, + implementation=implementation, + mode=mode, + targets=targets, + ) + ) + finally: + processor.dataset.close() + return rows + + +def run_target_probe( + *, + args: argparse.Namespace, + input_dir: Path, + run_root: Path, + target_sets: TargetSets, +) -> None: + target_detail_rows: list[dict[str, object]] = [] + for mode in args.modes: + target_detail_rows.extend( + matlab_target_rows( + run_root=run_root, + mode=mode, + skip_matlab=args.skip_matlab, + target_sets=target_sets, + ) + ) + target_detail_rows.extend( + python_target_rows( + args=args, + input_dir=input_dir, + run_root=run_root, + mode=mode, + target_sets=target_sets, + ) + ) + + target_summary_rows = aggregate_target_rows(target_detail_rows) + write_csv(run_root / "python_target_probe_detail.csv", target_detail_rows) + write_csv(run_root / "python_target_probe_summary.csv", target_summary_rows) + print(f"summary={run_root / 'python_target_probe_summary.csv'}") + print_table( + "MATLAB/Python vs implementation-specific fixed initial Python target", + target_summary_rows, + TARGET_TABLE_COLUMNS, + column_labels=TARGET_TABLE_LABELS, + notes=TARGET_TABLE_NOTES, + ) + + +def run_matlab_vs_python_probe( + *, + args: argparse.Namespace, + input_dir: Path, + run_root: Path, +) -> None: + detail_rows: list[dict[str, object]] = [] + for mode in args.modes: + matlab_dir = run_root / "matlab" / f"mode_{mode}" + if not matlab_dir.exists(): + raise FileNotFoundError(f"Missing MATLAB outputs for mode {mode}: {matlab_dir}") + for implementation in PYTHON_IMPLEMENTATIONS: + python_dir, processor = run_python_processor( + input_dir=input_dir, + output_root=run_root, + implementation=implementation, + mode=mode, + iterations=args.iterations, + seed=args.seed, + ) + processor.dataset.close() + detail_rows.extend(compare_dirs(matlab_dir, python_dir, mode, implementation)) + + summary_rows = aggregate(detail_rows) + write_csv(run_root / "matlab_shine_comparison_detail.csv", detail_rows) + write_csv(run_root / "matlab_shine_comparison_summary.csv", summary_rows) + print(f"summary={run_root / 'matlab_shine_comparison_summary.csv'}") + print_table( + "Python outputs relative to MATLAB SHINE", + summary_rows, + OUTPUT_COMPARISON_COLUMNS, + ) + + +def main() -> None: + args = parse_args() + run_root, input_dir, runner = prepare_run(args) + + target_sets = None + if args.use_python_targets: + print(ANSI_BOLD_RED + "\n".join(f"WARNING: {line}" for line in TARGET_WARNING_LINES) + ANSI_RESET) + target_sets = compute_target_sets( + input_dir=input_dir, + run_root=run_root, + seed=args.seed, + ) + if args.prepare_only: + return + if not args.skip_matlab: + run_matlab(args.matlab_bin, runner) + + input_gray_rows = compare_input_grayscale( + input_dir, + run_root, + seed=args.seed, + full_tracking=args.full_tracking, + ) + print_input_grayscale_results(input_gray_rows, run_root) + + if args.use_python_targets: + assert target_sets is not None + run_target_probe( + args=args, + input_dir=input_dir, + run_root=run_root, + target_sets=target_sets, + ) + else: + run_matlab_vs_python_probe(args=args, input_dir=input_dir, run_root=run_root) + + if not args.full_tracking: + cleanup_tracking_artifacts(run_root) + + +if __name__ == "__main__": + main() diff --git a/tests/tools/run_matlab_shine_comparison.sh b/tests/tools/run_matlab_shine_comparison.sh new file mode 100755 index 0000000..a86c32f --- /dev/null +++ b/tests/tools/run_matlab_shine_comparison.sh @@ -0,0 +1,329 @@ +#!/usr/bin/env bash +set -euo pipefail + +# Compare SHINIER (Python) against the original MATLAB SHINE toolbox. +# Requires MATLAB and the SHINE toolbox: +# http://www.mapageweb.umontreal.ca/gosselif/SHINE/ +# +# Centralized paths and run settings. Override any of these from the shell, e.g. +# RUN_ROOT=tmp/matlab_shine_comparison/my_run LIMIT=8 bash tests/tools/run_matlab_shine_comparison.sh +# If MATLAB/Python/SHINE are not found, the script asks for the needed paths. +SCRIPT_DIR="$(cd "$(dirname "${BASH_SOURCE[0]}")" && pwd)" +REPO_ROOT="$(cd "$SCRIPT_DIR/../.." && pwd)" + +PYTHON_BIN="${PYTHON_BIN:-python3}" +MATLAB_BIN="${MATLAB_BIN:-}" +SHINE_DIR="${SHINE_DIR:-}" + +OUTPUT_DIR="${OUTPUT_DIR:-$REPO_ROOT/tmp/matlab_shine_comparison}" +RUN_NAME="${RUN_NAME:-full_comparison}" +RUN_ROOT="${RUN_ROOT:-$OUTPUT_DIR/$RUN_NAME}" +MATLAB_VS_PYTHON_ROOT="$RUN_ROOT/matlab_vs_python" +PYTHON_TARGET_ROOT="$RUN_ROOT/python_targets" + +MODES="${MODES:-1 2 3 4 5 6 7 8}" +LIMIT="${LIMIT:-8}" +ITERATIONS="${ITERATIONS:-5}" +SEED="${SEED:-42}" +FULL_TRACKING="${FULL_TRACKING:-0}" + +PY_SCRIPT="$SCRIPT_DIR/matlab_shine_comparison.py" +REQUIRED_SHINE_FILES=(lumMatch.m histMatch.m sfMatch.m specMatch.m) +REQUIRED_SHINE_FUNCTIONS=(lumMatch histMatch sfMatch specMatch) + +read -r -a MODE_ARGS <<< "$MODES" + +BOLD="$(printf '\033[1m')" +RESET="$(printf '\033[0m')" + +is_executable() { + local candidate="$1" + command -v "$candidate" >/dev/null 2>&1 || [[ -x "$candidate" ]] +} + +detect_matlab_bin() { + local requested="$1" + local candidate + local found="" + + if [[ -n "$requested" && "$requested" != "matlab" && -x "$requested" ]]; then + printf '%s' "$requested" + return + fi + + if command -v matlab >/dev/null 2>&1; then + command -v matlab + return + fi + + for candidate in \ + /Applications/MATLAB_R*.app/bin/matlab \ + /Applications/MATLAB.app/bin/matlab \ + /usr/local/MATLAB/R*/bin/matlab \ + /opt/MATLAB/R*/bin/matlab \ + /usr/local/bin/matlab \ + /usr/bin/matlab \ + "/c/Program Files/MATLAB"/R*/bin/matlab.exe \ + "/mnt/c/Program Files/MATLAB"/R*/bin/matlab.exe + do + if [[ -x "$candidate" ]]; then + found="$candidate" + fi + done + + printf '%s' "$found" +} + +print_path_help() { + local label="$1" + + case "$label" in + Python) + printf '\n%sPython executable path%s\n' "$BOLD" "$RESET" >&2 + printf ' What to provide: the Python executable used to run this diagnostic.\n' >&2 + printf ' Examples:\n' >&2 + printf ' - python3\n' >&2 + printf ' - .venv/bin/python\n' >&2 + printf ' - /Users/you/environments/shinier_venv/bin/python\n' >&2 + ;; + MATLAB) + printf '\n%sMATLAB executable path%s\n' "$BOLD" "$RESET" >&2 + printf ' What to provide: the MATLAB executable, not the MATLAB.app folder.\n' >&2 + printf ' Examples:\n' >&2 + printf ' - /Applications/MATLAB_R2025a.app/bin/matlab\n' >&2 + printf ' - /usr/local/MATLAB/R2025a/bin/matlab\n' >&2 + printf ' - C:\\Program Files\\MATLAB\\R2025a\\bin\\matlab.exe\n' >&2 + ;; + SHINE) + printf '\n%sSHINE toolbox folder%s\n' "$BOLD" "$RESET" >&2 + printf ' What to provide: the folder containing MATLAB files such as:\n' >&2 + printf ' - lumMatch.m\n' >&2 + printf ' - histMatch.m\n' >&2 + printf ' - sfMatch.m\n' >&2 + printf ' - specMatch.m\n' >&2 + printf ' Example:\n' >&2 + printf ' - /Users/you/projects/shinetoolbox\n' >&2 + printf ' Download SHINE:\n' >&2 + printf ' - http://www.mapageweb.umontreal.ca/gosselif/SHINE/\n' >&2 + printf ' Press Enter only if MATLAB can already find those SHINE functions.\n' >&2 + ;; + esac +} + +prompt_required_path() { + local label="$1" + local current="$2" + local env_var="$3" + local value="$current" + + while ! is_executable "$value"; do + printf '%s not found: %s\n' "$label" "$value" >&2 + if [[ ! -t 0 ]]; then + print_path_help "$label" + printf 'Set %s=/path/to/executable or run this script interactively.\n' "$env_var" >&2 + exit 1 + fi + print_path_help "$label" + printf '\nEnter %s executable path: ' "$label" >&2 + read -r value + if [[ -z "$value" ]]; then + printf '%s path is required.\n' "$label" >&2 + exit 1 + fi + done + + printf '%s' "$value" +} + +prompt_optional_dir() { + local label="$1" + local current="$2" + local env_var="$3" + local value="$current" + + while [[ -n "$value" && ! -d "$value" ]]; do + printf '%s directory not found: %s\n' "$label" "$value" >&2 + if [[ ! -t 0 ]]; then + print_path_help "$label" + printf 'Set %s=/path/to/directory or leave it unset if MATLAB already has it on path.\n' "$env_var" >&2 + exit 1 + fi + print_path_help "$label" + printf '\nEnter %s directory path, or press Enter if MATLAB already has it on path: ' "$label" >&2 + read -r value + done + + if [[ -z "$value" && -t 0 && -z "${SHINE_DIR:-}" ]]; then + print_path_help "$label" + printf '\nEnter %s directory path, or press Enter if MATLAB already has it on path: ' "$label" >&2 + read -r value + while [[ -n "$value" && ! -d "$value" ]]; do + printf '%s directory not found: %s\n' "$label" "$value" >&2 + print_path_help "$label" + printf '\nEnter %s directory path, or press Enter if MATLAB already has it on path: ' "$label" >&2 + read -r value + done + fi + + printf '%s' "$value" +} + +validate_shine_dir() { + local missing=() + local file + + for file in "${REQUIRED_SHINE_FILES[@]}"; do + if [[ ! -f "$SHINE_DIR/$file" ]]; then + missing+=("$file") + fi + done + + if [[ "${#missing[@]}" -gt 0 ]]; then + printf '\n%sSHINE toolbox check failed%s\n' "$BOLD" "$RESET" >&2 + printf ' Folder checked:\n' >&2 + printf ' - %s\n' "$SHINE_DIR" >&2 + printf ' Missing required MATLAB files:\n' >&2 + for file in "${missing[@]}"; do + printf ' - %s\n' "$file" >&2 + done + printf ' Download SHINE:\n' >&2 + printf ' - http://www.mapageweb.umontreal.ca/gosselif/SHINE/\n' >&2 + printf ' Then rerun with SHINE_DIR=/path/to/shinetoolbox.\n' >&2 + exit 1 + fi +} + +validate_shine_matlab_path() { + local matlab_code + local function_list + + function_list="$(printf "'%s'," "${REQUIRED_SHINE_FUNCTIONS[@]}")" + function_list="${function_list%,}" + matlab_code="funcs={${function_list}}; missing={}; for k=1:numel(funcs), if exist(funcs{k}, 'file') ~= 2, missing{end+1}=funcs{k}; end; end; if ~isempty(missing), fprintf(2, '\\nSHINE toolbox check failed\\n'); fprintf(2, ' MATLAB cannot find required SHINE functions:\\n'); for k=1:numel(missing), fprintf(2, ' - %s\\n', missing{k}); end; fprintf(2, ' Provide SHINE_DIR=/path/to/shinetoolbox or download SHINE from:\\n'); fprintf(2, ' - http://www.mapageweb.umontreal.ca/gosselif/SHINE/\\n'); exit(1); end" + + if ! "$MATLAB_BIN" -batch "$matlab_code"; then + printf '\n%sSHINE toolbox check failed%s\n' "$BOLD" "$RESET" >&2 + printf ' MATLAB could not confirm that SHINE is on the MATLAB path.\n' >&2 + printf ' Provide SHINE_DIR=/path/to/shinetoolbox and rerun.\n' >&2 + exit 1 + fi +} + +validate_shine_setup() { + printf '\n[preflight] Checking SHINE toolbox availability\n' + if [[ -n "$SHINE_DIR" ]]; then + validate_shine_dir + else + validate_shine_matlab_path + fi +} + +run_matlab_runner() { + local runner="$1" + local matlab_runner="${runner//\'/\'\'}" + "$MATLAB_BIN" -batch "run('$matlab_runner')" +} + +cleanup_csv_only_root() { + local root="$1" + local path + + [[ -d "$root" ]] || return 0 + find "$root" -mindepth 1 -maxdepth 1 -type d -exec rm -rf {} + + while IFS= read -r -d '' path; do + case "$path" in + *.csv) ;; + *) rm -f "$path" ;; + esac + done < <(find "$root" -mindepth 1 -maxdepth 1 -type f -print0) +} + +cleanup_default_artifacts() { + if [[ "$FULL_TRACKING" == "1" || "$FULL_TRACKING" == "true" || "$FULL_TRACKING" == "yes" ]]; then + return 0 + fi + cleanup_csv_only_root "$MATLAB_VS_PYTHON_ROOT" + cleanup_csv_only_root "$PYTHON_TARGET_ROOT" +} + +print_config() { + printf '\nMATLAB SHINE vs SHINIER comparison\n' + printf '%-24s %s\n' 'repo_root' "$REPO_ROOT" + printf '%-24s %s\n' 'python_bin' "$PYTHON_BIN" + printf '%-24s %s\n' 'matlab_bin' "$MATLAB_BIN" + printf '%-24s %s\n' 'shine_dir' "${SHINE_DIR:-MATLAB path}" + printf '%-24s %s\n' 'run_root' "$RUN_ROOT" + printf '%-24s %s\n' 'modes' "$MODES" + printf '%-24s %s\n' 'limit' "$LIMIT" + printf '%-24s %s\n' 'iterations' "$ITERATIONS" + printf '%-24s %s\n' 'full_tracking' "$FULL_TRACKING" + printf '%-24s %s\n\n' 'seed' "$SEED" +} + +PYTHON_BIN="$(prompt_required_path Python "$PYTHON_BIN" PYTHON_BIN)" +MATLAB_BIN="$(detect_matlab_bin "$MATLAB_BIN")" +MATLAB_BIN="$(prompt_required_path MATLAB "$MATLAB_BIN" MATLAB_BIN)" +SHINE_DIR="$(prompt_optional_dir SHINE "$SHINE_DIR" SHINE_DIR)" + +COMMON_ARGS=( + --modes "${MODE_ARGS[@]}" + --limit "$LIMIT" + --iterations "$ITERATIONS" + --seed "$SEED" + --matlab-bin "$MATLAB_BIN" + --quiet-run-info +) + +if [[ -n "$SHINE_DIR" ]]; then + COMMON_ARGS+=(--shine-dir "$SHINE_DIR") +fi + +if [[ "$FULL_TRACKING" == "1" || "$FULL_TRACKING" == "true" || "$FULL_TRACKING" == "yes" ]]; then + COMMON_ARGS+=(--full-tracking) +fi + +cd "$REPO_ROOT" +mkdir -p "$RUN_ROOT" +trap cleanup_default_artifacts EXIT +print_config +validate_shine_setup + +printf '\n[1/6] Prepare MATLAB vs Python run\n' +"$PYTHON_BIN" "$PY_SCRIPT" "${COMMON_ARGS[@]}" \ + --prepare-only \ + --run-root "$MATLAB_VS_PYTHON_ROOT" + +printf '\n[2/6] Run MATLAB for MATLAB vs Python outputs\n' +run_matlab_runner "$MATLAB_VS_PYTHON_ROOT/run_matlab_shine.m" + +printf '\n[3/6] Compare MATLAB vs Python outputs\n' +"$PYTHON_BIN" "$PY_SCRIPT" "${COMMON_ARGS[@]}" \ + --skip-matlab \ + --run-root "$MATLAB_VS_PYTHON_ROOT" + +printf '\n[4/6] Prepare fixed Python target run\n' +"$PYTHON_BIN" "$PY_SCRIPT" "${COMMON_ARGS[@]}" \ + --use-python-targets \ + --prepare-only \ + --run-root "$PYTHON_TARGET_ROOT" + +printf '\n[5/6] Run MATLAB with fixed Python targets\n' +run_matlab_runner "$PYTHON_TARGET_ROOT/run_matlab_shine.m" + +printf '\n[6/6] Compare MATLAB/Python against fixed Python targets\n' +"$PYTHON_BIN" "$PY_SCRIPT" "${COMMON_ARGS[@]}" \ + --use-python-targets \ + --skip-matlab \ + --run-root "$PYTHON_TARGET_ROOT" + +printf '\nDone. Main outputs:\n' +printf '%-36s %s\n' 'MATLAB vs Python summary' "$MATLAB_VS_PYTHON_ROOT/matlab_shine_comparison_summary.csv" +printf '%-36s %s\n' 'MATLAB vs Python detail' "$MATLAB_VS_PYTHON_ROOT/matlab_shine_comparison_detail.csv" +printf '%-36s %s\n' 'Python target summary' "$PYTHON_TARGET_ROOT/python_target_probe_summary.csv" +printf '%-36s %s\n' 'Python target detail' "$PYTHON_TARGET_ROOT/python_target_probe_detail.csv" +printf '%-36s %s\n' 'Input grayscale comparison' "$MATLAB_VS_PYTHON_ROOT/input_grayscale_comparison.csv" +if [[ "$FULL_TRACKING" == "1" || "$FULL_TRACKING" == "true" || "$FULL_TRACKING" == "yes" ]]; then + printf '%-36s %s\n' 'Full tracking artifacts' "$RUN_ROOT" +else + printf '\nIntermediate images/scripts were removed. Set FULL_TRACKING=1 to keep them.\n' +fi diff --git a/tests/unit_tests/OptionalCconvolve_test.py b/tests/unit_tests/OptionalCconvolve_test.py new file mode 100644 index 0000000..6c87319 --- /dev/null +++ b/tests/unit_tests/OptionalCconvolve_test.py @@ -0,0 +1,61 @@ +import subprocess +import sys +import textwrap + +import pytest + + +pytestmark = pytest.mark.unit_tests + + +def test_import_falls_back_when_optional_cconvolve_fails(): + code = textwrap.dedent( + """ + import importlib.abc + import importlib.machinery + import sys + import warnings + + class BrokenCconvolveFinder(importlib.abc.MetaPathFinder): + def find_spec(self, fullname, path=None, target=None): + if fullname == "shinier._cconvolve": + return importlib.machinery.ModuleSpec(fullname, BrokenCconvolveLoader()) + return None + + class BrokenCconvolveLoader(importlib.abc.Loader): + def create_module(self, spec): + return None + + def exec_module(self, module): + raise ImportError("simulated NumPy ABI mismatch") + + sys.meta_path.insert(0, BrokenCconvolveFinder()) + + with warnings.catch_warnings(record=True) as caught: + warnings.simplefilter("always") + import shinier + + assert shinier._HAS_CYTHON is False + assert shinier.convolve2d_direct is None + assert shinier.convolve2d_separable is None + assert any( + "optional compiled convolution extension" in str(warning.message) + and "slower NumPy fallback" in str(warning.message) + and "Python:" in str(warning.message) + and "NumPy:" in str(warning.message) + and "reinstalling SHINIER" in str(warning.message) + and "simulated NumPy ABI mismatch" in str(warning.message) + for warning in caught + ) + """ + ) + + result = subprocess.run( + [sys.executable, "-W", "always", "-c", code], + stdout=subprocess.PIPE, + stderr=subprocess.PIPE, + text=True, + check=False, + ) + + assert result.returncode == 0, result.stdout + result.stderr