diff --git a/.github/workflows/tests.yml b/.github/workflows/tests.yml index 4908b96..30ec5f0 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 ] + branches: [ main, dev_1.4, dev_1.7, dev_v0.2.0, dev_v0.2.1] pull_request: - branches: [ main, dev_1.4, dev_1.7, dev_v0.2.0 ] + branches: [ main, dev_1.4, dev_1.7, dev_v0.2.0, dev_v0.2.1 ] jobs: tests: @@ -20,11 +20,11 @@ jobs: steps: # 1) Récupérer le code - name: Checkout repository - uses: actions/checkout@v4 + uses: actions/checkout@v5 # 2) Choisir la version de Python - name: Set up Python - uses: actions/setup-python@v5 + uses: actions/setup-python@v6 with: python-version: ${{ matrix.python-version }} diff --git a/README.md b/README.md index 4288cee..50bce5c 100644 --- a/README.md +++ b/README.md @@ -111,9 +111,8 @@ Below is an example of results obtained using mode 5 with joint histogram equali --- ## 🏛️ **Technical information** -See the accompanying the paper: [The SHINIER the Better: An Adaptation of the SHINE Toolbox on Python](documentation/papers/SHINIER/paper/paper.md) -And documentation: +Documentation: 1. [Package Overview](documentation/documentation.md#overview) 2. [Package Architecture](documentation/documentation.md#package-architecture) 3. [MATLAB vs Python Differences](documentation/documentation.md#matlab-vs-python-differences) @@ -133,7 +132,7 @@ And documentation: If you use **SHINIER**, please cite both of these articles: ### References -- 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* +- 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 --- diff --git a/documentation/demos.md b/documentation/demos.md index b7cdf85..3c63bfb 100644 --- a/documentation/demos.md +++ b/documentation/demos.md @@ -30,7 +30,7 @@ 6 = hist_match & spec_match 7 = sf_match & hist_match 8 = spec_match & hist_match -9 = only dithering +9 = standalone per-image transform (ie_methods or dithering) ``` --- @@ -332,21 +332,42 @@ opts = Options( --- -### Mode 9 – only dithering +### 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"`) ```python """ -Mode 9 (only dithering): applies noisy-bit dithering Allard & Faubert, 2008). +Mode 9 / ie_methods: applies a histogram equalization algorithm to each image + independently. Available algorithms: classic_he, tidhe, rdfhe, nfldice, betce, sfcef. -Example use case: dithering will be applied with the default noisy-bit method - (Allard & Faubert, 2008), while leaving the original image luminance and - spectrum unchanged. +Example use case: TIDHE is applied to each image independently. No inter-image + target is computed; each image is enhanced based on its own histogram. +""" +opts = Options( + input_folder=INPUT_FOLDER, + output_folder=OUTPUT_FOLDER, + mode = 9, + standalone_op = "ie_methods", + ie_methods = "tidhe", # classic_he | tidhe | rdfhe | nfldice | betce | sfcef +) +``` + +#### 9b – Dithering only (`standalone_op="dithering"`) + +```python +""" +Mode 9 / dithering: applies noisy-bit dithering (Allard & Faubert, 2008), + leaving the original image luminance and spectrum unchanged. """ opts = Options( input_folder=INPUT_FOLDER, output_folder=OUTPUT_FOLDER, mode = 9, - dithering = 1 + standalone_op = "dithering", + dithering = 1 # 1 = noisy-bit, 2 = Floyd-Steinberg ) ``` diff --git a/documentation/documentation.md b/documentation/documentation.md index 36f3849..60ad522 100644 --- a/documentation/documentation.md +++ b/documentation/documentation.md @@ -68,14 +68,7 @@ shinier/src ### Processing Flow -```mermaid -graph TD - A[Options] --> B[ImageDataset] - B --> C[ImageProcessor] - C --> D[Mode-based Processing] - D --> E[Final Conversion] - E --> F[Output Images] -``` +![SHINIER processing flow](figures/processing_flow.png) --- @@ -192,25 +185,48 @@ The choice between N and N-1 divisors depends on the **statistical context**: ### 4. **RGB to Grayscale Conversion** -#### MATLAB `rgb2gray()` - Legacy Standard +#### MATLAB `rgb2gray()` / NTSC-YIQ Intensity + +MATLAB-compatible grayscale conversion follows the Y channel of the historical +NTSC/YIQ transform. Although this conversion is often written with rounded +Rec.601-style coefficients, + ```matlab % Uses Rec.ITU-R BT.601-7 (SD monitors) -Y' = 0.299 * R + 0.587 * G + 0.114 * B +Y = 0.298936021293775 * R + 0.587043074451121 * G + 0.114020904255103 * B +``` + +In SHINIER, this MATLAB/NTSC-compatible intensity image can be obtained with: + +```python +gray = rgb2gray(image, weighting_standard="rec601", matlab_601=True) ``` -- **Rec.ITU-R BT.601-7**: Standard-Definition television (1982) -- **Legacy standard** optimized for CRT monitors + +or, equivalently: + +```python +gray = rgb2ntsc_intensity(image) +``` + +This path is mainly useful for MATLAB compatibility and NTSC/YIQ intensity +workflows. SHINIER's default modern grayscale preprocessing instead uses the CIE +`xyY` luminance channel after RGB linearization. When `legacy_mode=True` is set +in `Options`, SHINIER's color preprocessing automatically switches to this +MATLAB-compatible path (`rgb2gray(..., weighting_standard="rec601", matlab_601=True)`, +equivalent to `rgb2ntsc_intensity`) for grayscale images. #### SHINIER `rgb2gray()` - Modern Standards Support ```python -def rgb2gray(image, recommendation='rec709'): +def rgb2gray(image, weighting_standard='rec709', matlab_601=False): """RGB to grayscale conversion with multiple luma coefficient standards""" rgb_luma_coefficients = { 'equal': [0.333, 0.333, 0.333], # Equal weighting (not perceptually accurate) - 'rec601': [0.298936021293775, 0.587043074451121, 0.114020904255103], # SD legacy - 'rec709': [0.2126, 0.7152, 0.0722], # HD standard (recommended) - 'rec2020': [0.2627, 0.6780, 0.0593] # UHD standard + 'rec601': [0.222004309998231, 0.706654765925283, 0.0713409240764864], + 'rec709': [0.21263900587151, 0.715168678767756, 0.0721923153607337], + 'rec2020': [0.262700212011267, 0.677998071518871, 0.059301716469862], } - weights = rgb_luma_coefficients[recommendation] + matlab_rgb2gray_weights = [0.298936021293775, 0.587043074451121, 0.114020904255103] + weights = matlab_rgb2gray_weights if weighting_standard == 'rec601' and matlab_601 else rgb_luma_coefficients[weighting_standard] return np.dot(image.astype(np.float64), weights) ``` @@ -233,6 +249,12 @@ Different luma coefficients are optimized for different **display technologies** **WARNING AND REMINDER**: Ajusting for the luminance transfer functions implemented in image-capturing devices and the precise calibration of display monitors are essential for accurate visual stimuli presentation. +### 5. **Convolution — FMA and Unit in the Last Place (ULP)** + +MATLAB's `filter2` delegates to Intel MKL / Apple Accelerate, which uses **Fused Multiply-Add (FMA)** instructions. FMA computes `a × b + c` with a **single rounding step** instead of two, producing intermediate float64 values that can differ by ±1–2 ULP from NumPy's vectorised equivalent. + +At half-integer boundaries this shifts ~0.1 % of pixels by ±1 gray level. There is **no solution** in pure Python: the exact accumulation order inside MKL's `filter2` is undocumented and platform-dependent. + --- @@ -325,11 +347,26 @@ mode = 8 # spec_match → hist_match - Composite modes use temporary floating-point precision - Reduces rounding errors in multi-step calculations -### Mode 9: Dithering Only +### Mode 9: Standalone Per-Image Transform ```python -mode = 9 # only dithering +mode = 9 # standalone_op = "ie_methods" or "dithering" + +standalone_op = "ie_methods" # Image enhancement (default) +ie_methods = "tidhe" # Available algorithms: "tidhe", "rdfhe", "nfldice", "betce", "sfcef", "classic_he" + +standalone_op = "dithering" +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 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 @@ -440,6 +477,10 @@ class Options: fft_padding_mode: Literal[0, 1, 2, 3] = 0 fft_padding_value: Union[int, Literal[300]] = 300 + # --- Mode 9 --- + standalone_op: Literal["dithering", "ie_methods"] = "ie_methods" + ie_methods: Literal["classic_he", "tidhe", "rdfhe", "nfldice", "betce", "sfcef"] = "tidhe" + # --- Misc --- verbose: Literal[-1, 0, 1, 2, 3] = 0 ``` @@ -492,7 +533,7 @@ class ImageProcessor: verbose: Literal[-1, 0, 1, 2, 3] = 0 # --- Private attributes --- - _backward_conversion_type: str = PrivateAttr(default=None) + _backward_color_conversion: str = PrivateAttr(default=None) _color_space: Literal['uvw01', 'xyY'] = PrivateAttr(default='xyY') _complete: bool = PrivateAttr(default=False) _dataset_map: dict = PrivateAttr(default_factory=dict) @@ -687,6 +728,75 @@ mask_from_gui = masker.interactive_mask(image) 2. Quantize with rounding 3. Preserve overall image structure +### 5. Classic Global Histogram Equalization (Classic HE) + +**Algorithm:** +1. Compute the intensity histogram of the image +2. Compute the cumulative distribution function (CDF) +3. Map each intensity level *y* to `round(255 × CDF(y))` +4. Apply the look-up table to every pixel and cast to uint8 + +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) +**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:** +1. Find two partitioning levels where the cumulative histogram reaches ~1/3 and ~2/3 of total pixels (Eqs. 1-2) +2. Split the histogram into three equal-mass sub-bands: lower, middle, upper +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) +**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:** +1. Compute a fuzzy image histogram using the reference fuzziness parameter `p=10` +2. Find three recursive dualistic partitioning levels at ~25%, ~50%, and ~75% of fuzzy histogram mass +3. Split the fuzzy histogram into four sub-histograms +4. Equalize each fuzzy sub-histogram independently + +Implemented by `shinier.utils.rdfhe_gray`. + +### 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. + +**Algorithm:** +1. Fuzzify each gray level with the nonlinear (logistic) fuzzifier `F_X(I_o) = 1 / (1 + B^(−E_l·((I_o − P_l)/(L−1))))`, producing a membership value in (0, 1) (Eqs. 1-2) +2. Defuzzify the membership with the linear defuzzifier `D_L(I_f) = I_f × (L − 1)`, rescaling back to the gray-level range (Eqs. 3-4) +3. Apply the resulting monotonic look-up table to every pixel and cast to uint8 + +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) +**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. + +**Algorithm:** +1. Compute the entropy curve `ET(i) = -p_y(i) log2(p_y(i))` from the image intensity probabilities (Eq. 1) +2. Compute the partitioning level `pl` as the weighted arithmetic mean of gray levels using `ET(i)` as weights (Eq. 2) +3. Split `ET` into lower and upper sub-entropy curves `ET_l` and `ET_u` (Eqs. 3-4) +4. Equalize each sub-entropy curve independently and apply the resulting look-up table to every pixel (Eq. 5) + +Implemented by `shinier.utils.betce_gray`. + +### 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. + +**Algorithm:** +1. Set `phi=0.5`, `x=1`, and `t=0.5` for low-light images (`t=-1` is reported for low-contrast images) +2. Compute the coefficient bounds `a1`, `a2`, and `a3` of `G_S(phi)` (Eqs. 1-3) +3. Compute `c1` and `c2` from the fixed linear-combination weights `d1=d3=1/8`, `d2=d4=d6=0`, and `d5=1` (Eqs. 4-5) +4. Build the proposed 3x3 filter with `c1` around the center and `c2` at the center +5. Convolve the image with the 3x3 filter and cast to uint8 + +Implemented by `shinier.utils.sfcef_gray`. + --- @@ -727,18 +837,24 @@ class ImageListIO: - `Options`: Parameter validation - `ImageListIO`: Image loading/saving - `ImageDataset`: Collection management -- `ImageProcessor`: Image processing -- `Converter`: Luminance preservation and minimal chroma distortion. -- `GamutControl`: Chroma-loss minimization. +- `ImageProcessor`: Image processing pipeline and CLI +- `Converter`: Luminance preservation and minimal chroma distortion +- `GamutControl`: Chroma-loss minimization +- `Utils`: Utility functions (`rescale_images255`, histogram helpers, etc.) **Test Images:** - Noise-generated images for testing - Binary masks for figure-ground separation - Reference images for validation -### MATLAB Validation +### Validation Tests + +`ImageProcessor_validation_test.py` runs a combinatorial sweep over all `Options` parameters in three coverage modes (`sampled`, `pruned`, `exhaustive`). +It distinguishes **hard failures** (unexpected exceptions, broken internal invariants, SSIM final regression, RMSE more than doubled) from **soft failures** (minor RMSE regression expected in composite modes, SSIM sub-iteration rollback artifacts). +See `tests/README.md` for a detailed breakdown. -**TODO** +`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. --- diff --git a/documentation/figures/processing_flow.png b/documentation/figures/processing_flow.png new file mode 100644 index 0000000..344f5ff Binary files /dev/null and b/documentation/figures/processing_flow.png differ diff --git a/documentation/papers/SHINIER/paper/paper.bib b/documentation/papers/SHINIER/paper/paper.bib deleted file mode 100644 index 782112a..0000000 --- a/documentation/papers/SHINIER/paper/paper.bib +++ /dev/null @@ -1,62 +0,0 @@ -@article{allard2008noisy, - title={The noisy-bit method for digital displays: Converting a 256 luminance resolution into a continuous resolution}, - author={Allard, R{\'e}my and Faubert, Jocelyn}, - journal={Behavior Research Methods}, - volume={40}, - number={3}, - pages={735--743}, - year={2008}, - publisher={Springer}, - doi={10.3758/BRM.40.3.735} -} -@article{allen2022massive, - title={A massive 7T fMRI dataset to bridge cognitive neuroscience and artificial intelligence}, - author={Allen, Emily J and St-Yves, Ghislain and Wu, Yihan and Breedlove, Jesse L and Prince, Jacob S and Dowdle, Logan T and Nau, Matthias and Caron, Brad and Pestilli, Franco and Charest, Ian and others}, - journal={Nature neuroscience}, - volume={25}, - number={1}, - pages={116--126}, - year={2022}, - publisher={Nature Publishing Group US New York}, - doi={10.1038/s41593-021-00962-x} -} -@article{coltuc2006exact, - title={Exact histogram specification}, - author={Coltuc, Dinu and Bolon, Philippe and Chassery, J-M}, - journal={IEEE Transactions on Image processing}, - volume={15}, - number={5}, - pages={1143--1152}, - year={2006}, - publisher={IEEE}, - doi={10.1109/TIP.2005.864170} -} -@article{floyd1976adaptive, - title={An Adaptive Algorithm for Spatial Grayscale}, - author={Floyd, R. W. and Steinberg, L.}, - journal={Proceedings of the Society of Information Display}, - volume={17}, - number={2}, - pages={75--77}, - year={1976} -} -@article{srinath2017python, - title={Python--the fastest growing programming language}, - author={Srinath, KR}, - journal={International Research Journal of Engineering and Technology}, - volume={4}, - number={12}, - pages={354--357}, - year={2017} -} -@article{willenbockel2010controlling, - title={Controlling low-level image properties: the SHINE toolbox}, - author={Willenbockel, Verena and Sadr, Javid and Fiset, Daniel and Horne, Greg O and Gosselin, Fr{\'e}d{\'e}ric and Tanaka, James W}, - journal={Behavior research methods}, - volume={42}, - number={3}, - pages={671--684}, - year={2010}, - publisher={Springer}, - doi={10.3758/BRM.42.3.671} -} \ No newline at end of file diff --git a/documentation/papers/SHINIER/paper/paper.md b/documentation/papers/SHINIER/paper/paper.md deleted file mode 100644 index 791dc8f..0000000 --- a/documentation/papers/SHINIER/paper/paper.md +++ /dev/null @@ -1,93 +0,0 @@ ---- -title: 'The SHINIER the Better: An Adaptation of the SHINE Toolbox on Python' -tags: - - Python - - visual perception - - low-level image properties - - luminance - - histogram matching - - spatial frequency - - Fourier spectra -authors: - - name: Mathias Salvas-Hébert - orcid: 0009-0000-9707-7298 - corresponding: true - equal-contrib: true - affiliation: 1 - - name: Nicolas Dupuis-Roy - orcid: 0000-0001-9261-0583 - equal-contrib: true - affiliation: 2 - - name: Catherine Landry - orcid: 0000-0001-6748-1417 - affiliation: 1 - - name: Ian Charest - orcid: 0000-0002-3939-3003 - affiliation: 1 - - name: Frédéric Gosselin - orcid: 0000-0002-3797-4744 - affiliation: 1 -affiliations: - - name: Département de Psychologie, Université de Montréal, CP 6128, succ. Centre-ville, Montréal, QC, H3C 3J7, CANADA - index: 1 - - name: Elephant Scientific Consulting, Canada - index: 2 -date: 14 November 2025 -bibliography: paper.bib ---- - -# Summary - -SHINIER (Spectrum, Histogram, and Intensity Normalization, Equalization, and Refinement), written in Python, is an open-source package that replicates and extends the functionality of the popular SHINE toolbox [@willenbockel2010controlling], written in MATLAB. Like SHINE, it includes functions for normalizing and scaling mean luminance and contrast, for specifying either the full Fourier amplitude spectrum or its rotational average, and for exact histogram specification. In addition, SHINIER supports color images, better memory management, implements image dithering algorithms for improving pixel depth, and offers improved exact histogram equalization methods, among other enhancements. - -# Statement of need - -When conducting experiments with humans, animals, or machines, the choice of stimuli is critical. We usually intend observers to rely on features that genuinely support recognition in real life. However, experimental image sets—necessarily small subsets of the virtually infinite possible scenes and viewing conditions—often contain accidental features that can be exploited instead. For example, in a dog–cat categorization task, observers might succeed not because they attend to diagnostic shape or texture cues, but because the dog images (often taken outdoors in bright sunlight) have luminance histograms with higher means and greater variance than cat images (typically photographed indoors under dim lighting). These luminance differences are artifacts of illumination, not reliable distinguishing properties of dogs versus cats in the real world. One way to avoid such confounds is to use artificially generated stimuli with fully controlled low-level properties. Another is to rely on very large naturalistic image databases, such as the Natural Scenes Dataset [@allen2022massive], where idiosyncratic correlations tend to average out. When working with natural images and such large-scale resources are unavailable—or impractical due to time constraints with human or animal participants—normalizing and adjusting low-level image properties, as permitted by the SHINIER package, becomes essential. - -# State of the Field - -For more than fifteen years, the SHINE (Spectrum, Histogram, and Intensity Normalization and Equalization) toolbox [@willenbockel2010controlling] has served as the primary solution for controlling low-level image properties. It has been cited more than 1,400 times according to Google Scholar—an average of about 100 citations per year—clearly indicating its popularity and utility in the field. The original SHINE toolbox was written in MATLAB for controlling low-level image properties for vision research. Since then, the entire scientific ecosystem has shifted: Python has become the leading programming language of the scientific community [@srinath2017python]; computing capabilities have greatly improved; and vision science practices have evolved. - -To the best of our knowledge, no alternative open-source Python package has achieved comparable adoption or scope in vision science. General-purpose image processing Python libraries such as scikit-image or OpenCV are not designed to support the comprehensive control of low-level image properties required in vision research and addressed by SHINE. In particular, they do not provide exact histogram specification across image sets, spectrum normalization with preserved phase information, or workflows oriented toward reproducible stimulus generation. The SHINIER Python package (Spectrum, Histogram, and Intensity Normalization, Equalization, and Refinement) fills this gap by reimplementing the full functionality of SHINE in Python, while extending it to support current computational and experimental demands in vision science, including color image processing, high-precision numerical pipelines, and large-scale dataset handling. - -# Software Design - -## Accessibility - -Implementing the SHINIER package in Python enables seamless integration with modern data-science libraries (e.g., pandas, scikit-learn, PyTorch) that are not natively available in MATLAB. However, this choice also introduces legacy numerical compatibility issues stemming from intrinsic differences between Python and MATLAB. For example, the two languages rely on different rounding conventions. Numerical optimizations further increase discrepancies between the MATLAB and Python implementations. SHINIER therefore adopts a flexible numerical design that supports both legacy compatibility—through the legacy_mode option—and updated Python-based numerical handling. - -To maximize accessibility while simplifying long-term maintainability across versions, SHINIER was deliberately designed as a lightweight package with minimal dependencies. This design choice, however, required additional effort to preserve low-level performance through optimized implementations (e.g., using Cython for fast convolution), and to integrate external algorithms, such as a dedicated color subpackage for color-space support. - -## Simplified Structure - -SHINE followed the dominant script-based programming paradigm at the time of its conception. In contrast, SHINIER adopts an object-oriented programming (OOP) architecture that reorganizes the image-processing workflow into three core classes: (1) Options, a user-facing class that defines the input images and associated image processing parameters; (2) ImageDataset, which manages image storage and associated buffers; and (3) ImageProcessor, which applies the processing pipeline in a modular fashion. This modular organization improves maintainability and facilitates extension. - -## Improved Numerical Precision and Scalability - -Representing the image dataset as an object instantiation of a class provides additional opportunities for optimization and scalability. The numerical precision of images is now tracked across processing stages and, when possible, maintained in 64-bit floating-point format until the final step, where images are converted back to the mandatory 8-bit unsigned integer format. This eliminates precision bottlenecks inherited from the original SHINE workflow. - -The use of large datasets whose memory requirements exceed available random-access memory (RAM) is now the norm in computer vision and increasingly common in vision science. To ensure SHINIER’s scalability, we introduced ImageListIO, a dedicated subclass that can handle large datasets. This class leverages fast solid-state drive (SSD) read-and-write operations to circumvent RAM limitations. Specifically, when the conserve_memory option is enabled, images are read from and written on disk on the fly rather than stored in RAM. Although this new approach incurs lower real-time performance, in practice, the overhead is negligible compared to the core image-processing operations of the package. This approach also provides opportunities for future development, such as video processing. - -# New Functionalities and Algorithms - -SHINIER extends the functionality of SHINE in several ways. The package now supports color images and provides different luminance-processing options. When the linear_luminance option is enabled, input images are assumed to have RGB values linearly related to luminance; color channels are therefore processed independently, matching the original SHINE behavior for legacy compatibility. This approach can, however, produce out-of-gamut or distorted colors. We therefore provide alternative processing options that resolve this issue under different assumptions. - -The exact histogram workflow, which enforces identical pixel-intensity distributions across images, was completely revised. Two additional exact histogram specification algorithms were implemented: the method proposed by @coltuc2006exact and a Gaussian-kernel variant. New specification options were also introduced to handle special cases, including images with no identical pixel values—now common with high numerical precision—and situations in which a single algorithm is insufficient. In such cases, a hybrid strategy applies Gaussian ordering first and adds a tiny amount of noise only if ties remain. - -Finally, the Floyd–Steinberg dithering algorithm [@floyd1976adaptive] and the noisy-bit dithering algorithm [@allard2008noisy] were added to the package; in both cases, the effective pixel depth of the processed images is increased by exploiting spatial pooling in human and animal visual systems. - -# Research Impact Statement - -SHINIER should appeal to the hundreds of SHINE users migrating—or planning to migrate—from MATLAB to Python. It should also appeal to new users in the broader scientific community, including vision science and artificial intelligence (AI), who seek to control low-level image properties for applications involving humans, animals, or machines. Support for large image datasets that exceed available RAM, as well as for color images, further extends its applicability. Extensive unit and validation testing ensures numerical accuracy and reproducibility. - -The package is openly accessible via “pip install shinier”, with its source code hosted on Charestlab’s GitHub (https://github.com/Charestlab/SHINIER). Comprehensive documentation is provided in the README, and a command-line interface (CLI) allows user-friendly stimulus generation. - -# AI Usage Disclosure - -Different Large Language Models (LLM)—mainly ChatGPT, Claude, and Gemini—were used to assist in the following tasks: implementing the official style guide for Python code (PEP8); assisting in improving the OOP architecture; writing boilerplate code for classes and functions; providing insights for optimization; writing docstrings and API documentation; writing unit and validation tests. All proposed ideas and code were thoroughly examined, tested, and validated by human developers before committing them into the repository. - -# Acknowledgements - -The authors would like to thank the original contributors of the SHINE toolbox. - -# References \ No newline at end of file diff --git a/documentation/readthedocs/api.md b/documentation/readthedocs/api.md index a747cd4..4f7d3ad 100644 --- a/documentation/readthedocs/api.md +++ b/documentation/readthedocs/api.md @@ -16,7 +16,7 @@ documentation remains in the Markdown files under `documentation/`. :show-inheritance: .. autoclass:: shinier.ImageProcessor - :members: get_results, process, lum_match, hist_match, sf_match, spec_match + :members: get_results, process, lum_match, hist_match, sf_match, spec_match, dithering, ie_methods, uint8_to_float255, float255_to_float01, float01_to_float255, print_log_results, :exclude-members: __init__, __new__ :show-inheritance: ``` @@ -25,7 +25,7 @@ documentation remains in the Markdown files under `documentation/`. ```{eval-rst} .. autoclass:: shinier.color.ColorConverter - :members: sRGB_to_linRGB, linRGB_to_sRGB, linRGB_to_xyz, xyz_to_linRGB, xyz_to_xyY, xyY_to_xyz, xyz_to_lab, lab_to_xyz, sRGB_to_xyz, xyz_to_sRGB, sRGB_to_lab, lab_to_sRGB, sRGB_to_xyY, xyY_to_sRGB + :members: sRGB_to_linRGB, linRGB_to_sRGB, linRGB_to_xyz, xyz_to_linRGB, xyz_to_xyY, xyY_to_xyz, xyz_to_lab, lab_to_xyz, sRGB_to_xyz, xyz_to_sRGB, sRGB_to_lab, lab_to_sRGB, sRGB_to_xyY, xyY_to_sRGB, apply_standard_config :exclude-members: __init__, __new__ :show-inheritance: @@ -41,6 +41,8 @@ documentation remains in the Markdown files under `documentation/`. .. autofunction:: shinier.color.rgb2gray +.. autofunction:: shinier.color.rgb2ntsc_intensity + .. autofunction:: shinier.color.gray2rgb ``` @@ -61,9 +63,27 @@ documentation remains in the Markdown files under `documentation/`. .. autofunction:: shinier.utils.im_power_spectrum_plot +.. autofunction:: shinier.utils.tidhe_hist_plot + .. autofunction:: shinier.utils.show_processing_overview ``` +### Utility Classes + +```{eval-rst} +.. autoclass:: shinier.utils.Bcolors + :members: + +.. autoclass:: shinier.utils.DiffusionMaps + :members: + +.. autoclass:: shinier.utils.StepSizeController + :members: + +.. autoclass:: shinier.utils.ImageStats + :members: +``` + ### StimulusMasker ```{eval-rst} @@ -84,9 +104,9 @@ documentation remains in the Markdown files under `documentation/`. ```{eval-rst} .. autofunction:: shinier.utils.print_shinier_header -.. autofunction:: shinier.utils.get_field_values_from_pydantic_model +.. autofunction:: shinier.utils.pydantic_field_choices -.. autofunction:: shinier.utils.generate_pydantic_key_value_dict +.. autofunction:: shinier.utils.pydantic_model_choices .. autofunction:: shinier.utils.sf_profile @@ -182,17 +202,68 @@ documentation remains in the Markdown files under `documentation/`. .. autofunction:: shinier.utils.imhist +.. autofunction:: shinier.utils.classic_he_gray + +.. autofunction:: shinier.utils.tidhe_gray + +.. autofunction:: shinier.utils.rdfhe_gray + +.. autofunction:: shinier.utils.nfldice_gray + +.. autofunction:: shinier.utils.betce_gray + +.. autofunction:: shinier.utils.sfcef_gray + +``` + +### Contrast Enhancement Metrics + +```{eval-rst} +.. autofunction:: shinier.utils.compute_ambe + +.. autofunction:: shinier.utils.compute_contrast_improvement + +.. autofunction:: shinier.utils.compute_image_entropy + +.. autofunction:: shinier.utils.compute_mssim + +.. autofunction:: shinier.utils.compute_psnr + +.. autofunction:: shinier.utils.compute_bp2bpsim + +``` + +### Histogram And Image Statistics + +```{eval-rst} .. autofunction:: shinier.utils.rounded_target_hist .. autofunction:: shinier.utils.compute_tvd_hist .. autofunction:: shinier.utils.avg_hist + +.. autofunction:: shinier.utils.imstats ``` ## CLI ```{eval-rst} -.. autofunction:: shinier.SHINIER.SHINIER_CLI +.. autofunction:: shinier.SHINIER.prompt + +.. autofunction:: shinier.SHINIER.options_display + +.. autofunction:: shinier.SHINIER.get_image_list + +.. autofunction:: shinier.SHINIER.main +``` + +## Base Helpers + +```{eval-rst} +.. autofunction:: shinier.base.image_list_validator + +.. autoclass:: shinier.base.InformativeBaseModel + :exclude-members: __init__, __new__ ``` ## Others @@ -200,23 +271,27 @@ documentation remains in the Markdown files under `documentation/`. ### Chroma Loss Metrics ```{eval-rst} -.. autoclass:: shinier.color.quantify_chroma_loss.ChromaMetrics +.. autoclass:: shinier.color.chroma_eval.ChromaMetrics + +.. autoclass:: shinier.color.chroma_eval.AggregateMetric + +.. autoclass:: shinier.color.chroma_eval.AggregateRow -.. autoclass:: shinier.color.quantify_chroma_loss.AggregateMetric +.. autoclass:: shinier.color.chroma_eval.ChromaInfoRetention -.. autoclass:: shinier.color.quantify_chroma_loss.AggregateRow +.. autoclass:: shinier.color.chroma_eval.ChromaInfoLossResult -.. autoclass:: shinier.color.quantify_chroma_loss.ChromaInfoRetention +.. autofunction:: shinier.color.chroma_eval.chroma_metrics_for_image -.. autoclass:: shinier.color.quantify_chroma_loss.ChromaInfoLossResult +.. autofunction:: shinier.color.chroma_eval.aggregate_chroma_metrics -.. autofunction:: shinier.color.quantify_chroma_loss.compute_chroma_metrics_for_image +.. autofunction:: shinier.color.chroma_eval.build_chroma_loss_report -.. autofunction:: shinier.color.quantify_chroma_loss.aggregate_chroma_metrics +.. autofunction:: shinier.color.chroma_eval.chroma_info_loss_bpp -.. autofunction:: shinier.color.quantify_chroma_loss.generate_chroma_loss_report +.. autofunction:: shinier.color.chroma_eval.mean_chroma_loss_pct_lab -.. autofunction:: shinier.color.quantify_chroma_loss.chroma_info_loss_bits_per_pixel_vs_y1 +.. autofunction:: shinier.color.chroma_eval.lab_chroma_retention ``` ### ImageListIO diff --git a/src/shinier/ImageListIO.py b/src/shinier/ImageListIO.py index d6c47fb..230f29d 100644 --- a/src/shinier/ImageListIO.py +++ b/src/shinier/ImageListIO.py @@ -523,7 +523,7 @@ def _to_gray(self, image: np.ndarray) -> np.ndarray: """Convert an RGB image to grayscale according to the configured mode.""" if image.ndim == 3: if self.as_gray > 0: - image = rgb2gray(image, conversion_type=RGB2GRAY_WEIGHTS['int2key'][self.as_gray]) + image = rgb2gray(image, weighting_standard=RGB2GRAY_WEIGHTS['int2key'][self.as_gray]) image = uint8_plus(image) return image diff --git a/src/shinier/ImageProcessor.py b/src/shinier/ImageProcessor.py index 5280483..5a15e63 100644 --- a/src/shinier/ImageProcessor.py +++ b/src/shinier/ImageProcessor.py @@ -18,7 +18,9 @@ rescale_images255, get_images_spectra, ssim_sens, spectrum_plot, imhist_plot, sf_plot, avg_hist, uint8_plus, float01_to_uint, uint_to_float01, noisy_bit_dithering, floyd_steinberg_dithering, exact_histogram, Bcolors, MatlabOperators, compute_rmse, get_radius_grid, rotational_avg, - has_duplicates, stretch, console_log, print_log, StepSizeController, image_spectrum, _crop_after_fft + has_duplicates, stretch, console_log, print_log, StepSizeController, image_spectrum, _crop_after_fft, + IMAGE_ENHANCEMENT_METHODS, + compute_ambe, compute_contrast_improvement, compute_image_entropy, compute_mssim, compute_psnr, compute_bp2bpsim, ) from shinier.color import ColorConverter, ColorTreatment, rgb2gray, gray2rgb, RGB2GRAY_WEIGHTS, RGB_STANDARD, GamutControl @@ -102,13 +104,13 @@ class ImageProcessor(InformativeBaseModel): verbose: Literal[-1, 0, 1, 2, 3] = 0 # --- Private attributes --- - _backward_conversion_type: str = PrivateAttr(default=None) + _backward_color_conversion: str = PrivateAttr(default=None) _color_space: Literal['uvw01', 'xyY'] = PrivateAttr(default='xyY') _complete: bool = PrivateAttr(default=False) _dataset_map: dict = PrivateAttr(default_factory=dict) _fct_name2process_name: dict = PrivateAttr(default_factory=dict) _final_buffer: Optional[ImageListIO] = PrivateAttr(default=None) - _forward_conversion_type: str = PrivateAttr(default=None) + _forward_color_conversion: str = PrivateAttr(default=None) _gamut_control: Optional[GamutControl] = PrivateAttr(default=None) _initial_buffer: Optional[ImageListIO] = PrivateAttr(default=None) _initial_targets: Optional[Dict[str, np.ndarray]] = PrivateAttr(default={}) @@ -134,8 +136,8 @@ class ImageProcessor(InformativeBaseModel): def post_init(self, __context: Any) -> None: """Run initialization logic after Pydantic validation and only once at instantiation.""" - self._forward_conversion_type = f"sRGB_to_{self._color_space}" if self._color_space is not None else None - self._backward_conversion_type= f"{self._color_space}_to_sRGB" if self._color_space is not None else None + self._forward_color_conversion = f"sRGB_to_{self._color_space}" if self._color_space is not None else None + self._backward_color_conversion = f"{self._color_space}_to_sRGB" if self._color_space is not None else None if self.options is None: self.options = getattr(self.dataset, "options", None) @@ -165,7 +167,7 @@ def post_init(self, __context: Any) -> None: 6: ["hist_match", "spec_match"], 7: ["sf_match", "hist_match"], 8: ["spec_match", "hist_match"], - 9: [None], + 9: [None] if self.options.standalone_op == "dithering" else ["ie_methods"], } self._fct_name2process_name = { @@ -173,6 +175,7 @@ def post_init(self, __context: Any) -> None: "hist_match": "histogram matching", "sf_match": "spatial frequency matching", "spec_match": "fourier spectrum matching", + "ie_methods": "image enhancement", None: "dithering", } @@ -180,6 +183,10 @@ def post_init(self, __context: Any) -> None: self._processing_steps = self._mode2processing_steps[self.options.mode] self._n_steps = len(self._processing_steps) self._sum_bool_masks = [None] * len(self.dataset.images) + + if self.seed is None and self.options.seed is not None: + self.seed = self.options.seed + if self.from_unit_test: return @@ -203,29 +210,39 @@ def post_init(self, __context: Any) -> None: def print_log_results(self): """ - Compiles and prints logs related to regular processing steps, - SSIM validation results, and function validation results. The logs are - formatted and stored in the specified output folder. + Write selected options, regular processing logs, SSIM validation + results, and function validation results to the output log file. """ - # Print regular logs if self.verbose >= 0: - logs = [] + def format_option(value): + if isinstance(value, np.ndarray): + return f"ndarray(shape={value.shape}, dtype={value.dtype})" + if isinstance(value, Path): + return str(value) + return repr(value) + + def format_channel(value): + return "all" if value is None else value + + logs = ['[Options]'] + logs += [f"{key}: {format_option(value)}" for key, value in dict(self.options).items()] + logs += [''] + if len(self.log) > 0: logs = logs + ['[Regular logs]'] + self.log + [''] # Print SSIM validation results if len(self.ssim_results): - _logs = [f"iter={res['iter']}; step={res['step']}; image={res['image']}; channel={res['channel']}; result={res['valid_result']}" for res in self.ssim_results] + _logs = [f"iter={res['iter']}; step={res['step']}; image={res['image']}; channel={format_channel(res['channel'])}; result={res['valid_result']}" for res in self.ssim_results] logs = logs + ['[SSIM validation results]'] + _logs + [''] # Print function validation results if len(self.validation): - _logs = [f"iter={res['iter']}; step={res['step']}; image={res['image']}; channel={res['channel']}; processing function={res['processing_function']}; result={res['valid_result']}; other={res['log_result']}" for res in self.validation] + _logs = [f"iter={res['iter']}; step={res['step']}; image={res['image']}; channel={format_channel(res['channel'])}; processing function={res['processing_function']}; result={res['valid_result']}; other={res['log_result']}" for res in self.validation] logs = logs + ['[Processing validation results]'] + _logs + [''] - if len(logs): - print_log(logs=logs, log_path=self.options.output_folder) + print_log(logs=logs, log_path=self.options.output_folder) def get_results(self): """Return the processed output in the same container style as the input. @@ -354,7 +371,7 @@ def _compute_target_spectrum_from_image_path(self, image_path: Path) -> np.ndarr linear_luminance=self.options.linear_luminance, as_gray=self.options.as_gray, output_other=target_buffer_other, - conversion_type=self._forward_conversion_type, + color_conversion=self._forward_color_conversion, desaturate_chroma_on_low_luminance=True, legacy_mode=self.options.legacy_mode, verbose=False, @@ -396,7 +413,7 @@ def _compute_initial_target_histogram(self, n_bins: int = 256): 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: - raise ValueError(f"target_hist must have {self.dataset.buffer.n_channels} channels, ") + raise ValueError(f"target_hist must have {self.dataset.buffer.n_channels} channels, but has {target_hist.shape[-1]}.") self._target_hist = target_hist def _compute_target_hist_from_image_path(self, image_path: Path, n_bins: int = 256) -> np.ndarray: @@ -433,7 +450,7 @@ def _compute_target_hist_from_image_path(self, image_path: Path, n_bins: int = 2 linear_luminance=self.options.linear_luminance, as_gray=self.options.as_gray, output_other=target_buffer_other, - conversion_type=self._forward_conversion_type, + color_conversion=self._forward_color_conversion, desaturate_chroma_on_low_luminance=True, legacy_mode=self.options.legacy_mode, verbose=False, @@ -503,8 +520,8 @@ def _prepare_mask(image, mask=None): # Compute number of unmasked pixels per channel self._sum_bool_masks[idx] = [self.bool_masks[idx][..., ch].sum() for ch in range(self.bool_masks[idx].shape[2])] - def _validate_ssim(self, ssim: List[float]): - """Validate that SSIM progression is monotonic during optimization.""" + def _validate_ssim(self, ssim: List[float], tag: str = 'sub_iter'): + """Validate that SSIM did not regress. tag: 'sub_iter' (pre-rollback) or 'final' (output).""" out = np.array(ssim) if out.shape[0] > 1: for ch in range(out.shape[1]): @@ -514,15 +531,18 @@ def _validate_ssim(self, ssim: List[float]): 'step': self._step, 'image': self._processed_image, 'channel': self._processed_channel, - 'valid_result': is_strictly_increasing + 'tag': tag, + 'valid_result': is_strictly_increasing, } + if not is_strictly_increasing: + results['ssim_values'] = out[:, ch].tolist() self.ssim_results.append(results) if not is_strictly_increasing and self.verbose > 1: - res = f'{Bcolors.OKCYAN}SSIM optimization test for channel {ch}:{Bcolors.ENDC} {Bcolors.FAIL}FAIL{Bcolors.ENDC}' + res = f'{Bcolors.OKCYAN}SSIM optimization {tag} test for channel {ch}:{Bcolors.ENDC} {Bcolors.FAIL}FAIL{Bcolors.ENDC}' console_log(msg=res, indent_level=1, verbose=self.verbose > 2) - raise Exception(f"SSIM optimization non-monotonic for channel {ch}: {out[:, ch]}") + raise Exception(f"SSIM optimization {tag} non-monotonic for channel {ch}: {out[:, ch]}") - res = f'{Bcolors.OKCYAN}SSIM optimization test:{Bcolors.ENDC} {Bcolors.OKGREEN}PASS{Bcolors.ENDC}' + res = f'{Bcolors.OKCYAN}SSIM optimization {tag} test:{Bcolors.ENDC} {Bcolors.OKGREEN}PASS{Bcolors.ENDC}' console_log(msg=res, indent_level=1, verbose=self.verbose >= 3) def _validate(self, observed: List[float], expected: List[float], measures_str: list[str], rmse_tolerance: float = 1e-3): @@ -551,6 +571,38 @@ def _validate(self, observed: List[float], expected: List[float], measures_str: console_log(msg=results['log_result'], indent_level=indent_level, verbose=self.verbose==3) self.validation.append(results) + def _record_ie_metrics(self, before: np.ndarray, after: np.ndarray, op_label: str) -> None: + """Record image-enhancement metrics as a diagnostic validation entry.""" + before_ci = compute_contrast_improvement(before) + after_ci = compute_contrast_improvement(after) + before_entropy = compute_image_entropy(before) + after_entropy = compute_image_entropy(after) + warning = after_ci < before_ci and after_entropy < before_entropy + result = "WARN" if warning else "PASS" + metrics = ( + f"AMBE={compute_ambe(before, after):.6f}; " + f"MSSIM={compute_mssim(before, after):.6f}; " + f"PSNR={compute_psnr(before, after):.6f}; " + f"BP2BPSIM={compute_bp2bpsim(before, after):.6f}; " + f"CI={before_ci:.6f}->{after_ci:.6f} (delta={after_ci - before_ci:.6f}); " + f"Entropy={before_entropy:.6f}->{after_entropy:.6f} (delta={after_entropy - before_entropy:.6f})" + ) + note = ( + "WARN: CI and entropy both decreased; the image may already be bright/high-contrast, " + "or the transform may be reducing useful contrast." + if warning else + "PASS: diagnostic metrics did not trigger the CI+entropy decrease warning." + ) + self.validation.append({ + 'iter': self._iter_num, + 'step': self._step, + 'processing_function': "ie_methods", + 'image': self._processed_image, + 'channel': "all", + 'valid_result': result, + 'log_result': f"{op_label}: {metrics}; {note}", + }) + def process(self): """Run the full SHINIER pipeline on the current dataset. @@ -610,7 +662,7 @@ def process(self): linear_luminance=self.options.linear_luminance, as_gray=self.options.as_gray, output_other=self.dataset.buffer_other, - conversion_type=self._forward_conversion_type, + color_conversion=self._forward_color_conversion, desaturate_chroma_on_low_luminance=True, legacy_mode=self.options.legacy_mode, verbose=self.verbose>=2) @@ -711,7 +763,8 @@ def process(self): linear_luminance=self.options.linear_luminance, gamut_strategy=self.options.gamut_strategy, as_gray=self.options.as_gray, - conversion_type=self._backward_conversion_type, + color_conversion=self._backward_color_conversion, + legacy_mode=self.options.legacy_mode, verbose=self.verbose>=2) # Applies dithering or simply convert into uint8 if no dithering @@ -844,8 +897,8 @@ def compute_stats(im: np.ndarray, binary_mask: np.ndarray) -> Tuple[float, float # M = MatlabOperators.mean2(im[binary_mask]) if self.options.legacy_mode else np.mean(im[binary_mask]) # SD = MatlabOperators.mean2(im[binary_mask]) if self.options.legacy_mode else np.mean(im[binary_mask]) # else: - # convertion_type = RGB2GRAY_WEIGHTS['int2key'][self.options.rgb_weights] - # ch_weights = RGB2GRAY_WEIGHTS[conversion_type] + # weighting_standard = RGB2GRAY_WEIGHTS['int2key'][self.options.rgb_weights] + # ch_weights = RGB2GRAY_WEIGHTS[weighting_standard] # ch_means = np.array([np.mean(im[:, :, c][binary_mask[:, :, c]]) for c in range(3)]) # ch_stds = np.array([np.std(im[:, :, c][binary_mask[:, :, c]]) for c in range(3)]) # M = np.sum(ch_means * ch_weights) @@ -1196,7 +1249,7 @@ def hist_match(self): Y, _ = exact_histogram(image=X, binary_mask=self.bool_masks[idx], target_hist=self._target_hist, tie_strategy='none', n_bins=n_bins) else: Y, OA = exact_histogram(image=X, binary_mask=self.bool_masks[idx], target_hist=self._target_hist, tie_strategy=tie_strategy, n_bins=n_bins) - if hist_spec_names != 'noise' and (n_iter == 1 or (n_iter > 1 and self._sub_iter < n_iter - 1)): + if tie_strategy != 'noise' and (n_iter == 1 or (n_iter > 1 and self._sub_iter < n_iter - 1)): console_log(msg=f"Ordering accuracy per channel = {OA}", indent_level=1, color=Bcolors.OKBLUE, verbose=self.verbose == 3) # Compute Structural Similarity and gradient map (sens), along with max and min if self._sub_iter < n_iter - 1: @@ -1244,9 +1297,14 @@ def hist_match(self): self._sub_iter += 1 - # Test monotonic increase of ssim between first and last iteration - if self.options.hist_optim and len(all_ssim) >=2: - self._validate_ssim(ssim=[all_ssim[0], all_ssim[-1]]) + if self.options.hist_optim and len(all_ssim) >= 2: + # sub_iter: first vs last pre-rollback proposal — shows per-step trajectory. + self._validate_ssim(ssim=[all_ssim[0], all_ssim[-1]], tag='sub_iter') + if self.options.hist_optim and len(all_ssim) >= 1 and self.from_validation_test: + # final: first-iteration SSIM vs actual output Y — definitive correctness check. + # Guarded by from_validation_test: ssim_sens is O(W×H) and too costly in production. + _, final_ssim = ssim_sens(original_image, Y, data_range=n_bins - 1, use_sample_covariance=False, binary_mask=self.bool_masks[idx]) + self._validate_ssim(ssim=[all_ssim[0], final_ssim], tag='final') # Important Note: # - Must use Y as this is the one that matches the target histogram. @@ -1266,6 +1324,88 @@ def hist_match(self): buffer_collection.drange = (0, 255) + def ie_methods(self): + """Apply standalone image enhancement to each image independently. + + Dispatches based on ``options.ie_methods``: + + - ``'classic_he'``: Classic global histogram equalization. Maps each intensity + level *y* to ``round(255 × CDF(y))``, producing an approximately flat + output histogram. Maximizes contrast but can over-enhance noise on natural images. + Implemented by :func:`shinier.utils.classic_he_gray`. + + - ``'tidhe'``: Tripartite Image Decomposition-Based Histogram Equalization + (Rahman & Shimamura, *ICIC Express Letters* 20(3), 2026). Splits the + histogram into three equal-mass sub-bands, clips each via Intensity + Histogram Clipping (IHC), and equalizes them independently. + Implemented by :func:`shinier.utils.tidhe_gray`. + + - ``'rdfhe'``: Recursive Dualistic Fuzzy Histogram Equalization + (Rahman et al., IEEE QPAIN, 2026). Computes a fuzzy histogram, + partitions it recursively into four sub-histograms, and equalizes them + independently. Implemented by :func:`shinier.utils.rdfhe_gray`. + + - ``'nfldice'``: Nonlinear Fuzzification–Linear Defuzzification-Based ICE + (Rahman, *Trends in Electronics and Health Informatics*, LNNS 1034, 2025). + A fuzzy set-theoretic point operation: fuzzifies each gray level with a + nonlinear (logistic) fuzzifier, then defuzzifies linearly back to the + gray-level range. Implemented by :func:`shinier.utils.nfldice_gray`. + + - ``'betce'``: Bi-Entropy Curve Equalization + (Rahman, IEEE EICT, 2025). A state-of-the-art curve-based algorithm: + replaces the histogram with an entropy curve, partitions it into two + sub-curves, then equalizes each independently. + Implemented by :func:`shinier.utils.betce_gray`. + + - ``'sfcef'``: Sakaguchi-type Function-Based Cost-Effective Filtering + (Rahman et al., *Pattern Analysis and Applications*, 2025). A + state-of-the-art filtering-based algorithm: builds one 3x3 filter from + Sakaguchi/Gegenbauer coefficient bounds and convolves each image. + Implemented by :func:`shinier.utils.sfcef_gray`. + + In all cases each image is processed independently — no shared target is + computed across the dataset. The transform is applied per channel on the + prepared luminance buffer (or on all channels when + ``linear_luminance=True``). + + Results are stored back in ``dataset.buffer`` as ``float64`` in ``[0, 255]``. + The final output conversion casts these values to ``uint8`` without + dithering; mode 9 image enhancement does not accept a dithering option. + """ + buffer_collection = self.dataset.buffer + n_bins = 256 + op = self.options.ie_methods + ie_method = IMAGE_ENHANCEMENT_METHODS[op] + ie_fn = ie_method["fn"] + op_label = ie_method["label"] + + for idx, image in enumerate(buffer_collection): + self._processed_image = f'#{idx}' if self.dataset.images.src_paths[idx] is None else self.dataset.images.src_paths[idx] + console_log(msg=f"\nImage {self._processed_image}", indent_level=0, color=Bcolors.BOLD, verbose=self.verbose >= 2) + + X = im3D(image) + X_int = MatlabOperators.uint8(X) if self.options.legacy_mode \ + else np.clip(np.round(X), 0, 255).astype(np.uint8) + Y = X.copy() + for ch in range(X_int.shape[-1]): + Y[..., ch] = ie_fn(X_int[..., ch], legacy_mode=self.options.legacy_mode) + + self._record_ie_metrics(X_int, Y, op_label) + + final_hist = imhist(image=Y.squeeze() if Y.shape[-1] == 1 else Y, mask=self.bool_masks[idx], n_bins=n_bins, normalized=True) + equal_target = np.ones(final_hist.shape) / n_bins + rmse = compute_rmse(final_hist.flatten(), equal_target.flatten()) + console_log( + msg=f"Histogram flatness error ({op_label}): {rmse:.5f}", + indent_level=1, + color=Bcolors.OKBLUE, + verbose=self.verbose == 3, + ) + + buffer_collection[idx] = Y + + buffer_collection.drange = (0, 255) + def sf_match(self): """Match the rotationally averaged spatial-frequency content of each image (i.e., mean amplitude per spatial frequency). @@ -1349,9 +1489,10 @@ def sf_match(self): rmse = compute_rmse(t, o) self._validate(observed=[rmse], expected=[0], measures_str=['RMS error']) - # Soft-clip output values: As this transformation typically produces out-of-range values + # Soft-clip on last outer iteration: as this transformation typically produces out-of-range values output_image = np.stack(matched_image, axis=-1).squeeze() - if self._is_last_operation: + # Last iteration to include iterative modes with sf_match before the last operation (e.g., mode 5 or 7) + if self._iter_num == self.options.iterations - 1: mn, mx = output_image.min(), output_image.max() if mn < 0 or mx > 1: console_log( @@ -1442,9 +1583,10 @@ def spec_match(self): self._validate(observed=[rmse], expected=[0], measures_str=['RMS error']) - # Soft-clip output values: As this transformation typically produces out-of-range values + # Soft-clip on last outer iteration: as this transformation typically produces out-of-range values output_image = np.stack(matched_image, axis=-1).squeeze() - if self._is_last_operation: + # Last iteration to include iterative modes with spec_match before the last operation (e.g., mode 6 or 8) + if self._iter_num == self.options.iterations - 1: mn, mx = output_image.min(), output_image.max() if mn < 0 or mx > 1: console_log( diff --git a/src/shinier/Options.py b/src/shinier/Options.py index d995418..fecd25a 100644 --- a/src/shinier/Options.py +++ b/src/shinier/Options.py @@ -33,6 +33,7 @@ 'luminance': ['safe_lum_match', 'target_lum'], 'histogram': ['hist_optim', 'hist_specification', 'hist_iterations', 'target_hist'], 'fourier': ['rescaling', 'target_spectrum', 'fft_padding_mode', 'fft_padding_value'], + 'mode9': ['standalone_op', 'ie_methods'], 'misc': ['verbose'] } @@ -59,7 +60,9 @@ class Options(InformativeBaseModel): ``hist_optim``, ``hist_specification``, ``hist_iterations``, ``target_hist`` 8. FOURIER matching ``rescaling``, ``target_spectrum``, ``fft_padding_mode``, ``fft_padding_value`` - 9. Misc + 9. MODE 9 + ``standalone_op``, ``ie_methods`` + 10. Misc ``verbose`` Options @@ -115,7 +118,7 @@ class Options(InformativeBaseModel): - 6 = hist_match and spec_match. - 7 = sf_match and hist_match. - 8 = spec_match and hist_match. - - 9 = only dithering. + - 9 = standalone per-image transform — no inter-image target, controlled by ``standalone_op``: Related methods in :class:`shinier.ImageProcessor`: @@ -123,6 +126,7 @@ class Options(InformativeBaseModel): - Histogram matching: :meth:`shinier.ImageProcessor.hist_match` - Spatial-frequency matching: :meth:`shinier.ImageProcessor.sf_match` - Spectrum matching: :meth:`shinier.ImageProcessor.spec_match` + - Image enhancement: :meth:`shinier.ImageProcessor.ie_methods` - Full pipeline orchestration: :meth:`shinier.ImageProcessor.process` legacy_mode : Optional[bool] @@ -133,17 +137,23 @@ class Options(InformativeBaseModel): Important: legacy_mode affects more than the explicit option overrides listed below. It also enables MATLAB-compatibility behavior in several processing steps - (for example MATLAB-style rounding and grayscale conversion paths), so outputs may differ - even when two runs appear to share the same visible option values. + (for example MATLAB-style rounding and MATLAB-compatible RGB to grayscale + conversion), so outputs may differ even when two runs appear to share the + same visible option values. True reproduces the behavior of previous releases by setting: - ``conserve_memory`` = ``False`` - ``as_gray`` = ``1`` + - ``linear_luminance`` = ``False`` + - ``rec_standard`` = ``1`` (Rec.601) - ``dithering`` = ``0`` - ``hist_specification`` = ``1`` - ``safe_lum_match`` = ``False`` + In grayscale color preprocessing, legacy mode uses MATLAB-compatible + Rec.601 / NTSC-YIQ intensity weights, equivalent to MATLAB ``rgb2gray``. + False means no legacy settings are forced and all options follow their current defaults. seed : Optional[int] @@ -158,7 +168,7 @@ class Options(InformativeBaseModel): iterations : int [3] SHINIER MODE. - Number of iterations for composite modes. + Number of iterations for composite modes (5-8). Default is 5. For these modes, histogram specification and Fourier amplitude specification affect each other. @@ -167,6 +177,8 @@ class Options(InformativeBaseModel): This method of iterating was developed so that it recalculates the respective target at each iteration (i.e., no target hist/spectrum). + Silently forced to 1 outside composite modes (1-4, 9). + as_gray : bool [4] Grayscale / color. @@ -176,8 +188,11 @@ class Options(InformativeBaseModel): - True = Convert into grayscale images. - When ``linear_luminance`` is ``False``: - computes non-linear grayscale images by applying the perceptual - luma weights from the specified ``rec_standard``. + extracts the CIE xyY luminance channel after RGB + linearization, using the selected ``rec_standard`` for the + RGB-to-XYZ conversion. In ``legacy_mode``, SHINIER instead + uses MATLAB-compatible Rec.601 ``rgb2gray`` / NTSC-YIQ + intensity weights. - When ``linear_luminance`` is ``True``: computes linear grayscale images by averaging the RGB channels @@ -364,7 +379,9 @@ class Options(InformativeBaseModel): - Spatial dimensions must match the processed images. - ``'equal'``: - - Uses a flat histogram, i.e., histogram equalization. + - Flat histogram, i.e., performs a type of `histogram equalization`. + - All images are matched to the same flat histogram (``1/256`` per bin). + - Preserves inter-image consistency; uses the full ``exact_histogram`` pipeline. - ``None``: - Uses the average histogram of all input images. @@ -378,11 +395,11 @@ class Options(InformativeBaseModel): Default is 2. - 0 = no rescaling. - - 1 = Rescaling each image so that it stretches to [0, 1] (its own min -> 0, max -> 1). - - 2 = Rescaling absolute max/min (shared 0-1 range). - - 3 = Rescaling average max/min. + - 1 = Per-image stretch to [0, 255] (each image's own min → 0, max → 255). + - 2 = Dataset absolute min/max mapped to [0, 255] (shared range, no clipping). + - 3 = Dataset average min/max mapped to [0, 255] (shared range, outlier images are clipped). - Not allowed for modes 1 and 2. + Not used in modes 1 and 2 (silently reset to 0). target_spectrum : Optional[Union[np.ndarray, Path]] [8] FOURIER matching. @@ -430,8 +447,44 @@ class Options(InformativeBaseModel): If ``300``, the mean intensity of the current normalized image is used. Used only when ``fft_padding_mode=3``. + Silently reset to ``300`` unless ``fft_padding_mode=3``. + + standalone_op : Literal['dithering', 'ie_methods'] + [9] MODE 9. + + Selects the standalone transform applied in mode 9. + Default is ``'ie_methods'``. + + Only used when ``mode=9``. + + - ``'dithering'``: applies the dithering method selected by ``dithering`` (1 or 2) before + the final uint8 cast. Requires ``dithering != 0``. + - ``'ie_methods'``: applies image enhancement per image; the specific + algorithm is controlled by ``ie_methods``. + + ie_methods : Literal['classic_he', 'tidhe', 'rdfhe', 'nfldice', 'betce', 'sfcef'] + [9] MODE 9. + + Selects the image-enhancement algorithm applied when ``standalone_op='ie_methods'``. + Default is ``'tidhe'``. + + Only used when ``mode=9`` and ``standalone_op='ie_methods'``. + + - ``'classic_he'``: classic global histogram equalization. + See :func:`shinier.utils.classic_he_gray`. + - ``'tidhe'``: Tripartite Image Decomposition-Based Histogram Equalization. + See :func:`shinier.utils.tidhe_gray`. + - ``'rdfhe'``: Recursive Dualistic Fuzzy Histogram Equalization. + See :func:`shinier.utils.rdfhe_gray`. + - ``'nfldice'``: Nonlinear Fuzzification-Linear Defuzzification-Based ICE. + See :func:`shinier.utils.nfldice_gray`. + - ``'betce'``: Bi-Entropy Curve Equalization. + See :func:`shinier.utils.betce_gray`. + - ``'sfcef'``: Sakaguchi-type Function-Based Cost-Effective Filtering. + See :func:`shinier.utils.sfcef_gray`. + verbose : Literal[-1, 0, 1, 2, 3], optional - [9] Misc. + [10] Misc. Controls verbosity levels. Default is 0. @@ -490,6 +543,10 @@ class Options(InformativeBaseModel): fft_padding_mode: Literal[0, 1, 2, 3] = 0 fft_padding_value: Union[int, Literal[300]] = 300 + # --- Mode 9 --- + standalone_op: Literal["dithering", "ie_methods"] = "ie_methods" + ie_methods: Literal["classic_he", "tidhe", "rdfhe", "nfldice", "betce", "sfcef"] = "tidhe" + # --- Misc --- verbose: Literal[-1, 0, 1, 2, 3] = 0 @@ -512,27 +569,32 @@ def validate_existing_path(cls, v: Optional[Path]) -> Optional[Path]: @field_validator("target_hist") @classmethod def validate_target_hist(cls, v): - """Validate that target_hist is 'equal', an array of correct shape, or a valid image path.""" - if v is None or (isinstance(v, str) and v in ["equal", 'unit_test']): + """Validate that target_hist is 'equal', a numpy array of correct shape, or an image path.""" + if v is None: return v - if isinstance(v, (str, Path)): - v = Path(v).resolve() - if not v.exists(): - raise ValueError(f"target_hist image does not exist: {v}") - if not v.is_file(): - raise ValueError(f"target_hist path must point to a file: {v}") - if v.suffix.lower().lstrip(".") not in get_args(ACCEPTED_IMAGE_FORMATS): - raise ValueError( - f"target_hist image must use one of {get_args(ACCEPTED_IMAGE_FORMATS)}. " - f"Got: {v.suffix}" - ) + if isinstance(v, np.ndarray): + if v.ndim not in (1, 2): + raise ValueError("target_hist must be 1D (gray) or 2D (color).") + if v.ndim == 1 and v.size != 256: + raise ValueError("For grayscale, target_hist must have 256 bins.") return v - if not isinstance(v, np.ndarray): + if isinstance(v, (str, Path)) and str(v) in {"equal", "unit_test"}: + return str(v) + if isinstance(v, str): + raise ValueError(f"target_hist string value '{v}' is not recognised. " + f"Use 'equal', an image path, or a numpy array.") + if not isinstance(v, Path): raise TypeError("target_hist must be a numpy.ndarray, an image path, or 'equal'.") - if v.ndim not in (1, 2): - raise ValueError("target_hist must be 1D (gray) or 2D (color).") - if v.ndim == 1 and v.size != 256: - raise ValueError("For grayscale, target_hist must have 256 bins.") + v = v.resolve() + if not v.exists(): + raise ValueError(f"target_hist image does not exist: {v}") + if not v.is_file(): + raise ValueError(f"target_hist path must point to a file: {v}") + if v.suffix.lower().lstrip(".") not in get_args(ACCEPTED_IMAGE_FORMATS): + raise ValueError( + f"target_hist image must use one of {get_args(ACCEPTED_IMAGE_FORMATS)}. " + f"Got: {v.suffix}" + ) return v @field_validator("target_spectrum") @@ -601,9 +663,18 @@ def cross_checks(self) -> "Options": object.__setattr__(self, "rescaling", 0) console_log(msg=f"Rescaling not valid for luminance/histogram modes. rescaling -> 0", color=Bcolors.WARNING, verbose=self.verbose > 0) - # Mode 9: must have dithering != 0 → raise ValueError - if self.mode == 9 and self.dithering == 0: - raise ValueError("Mode 9 requires dithering 1 or 2 (not 0).") + # Mode 9 (dithering): must have dithering != 0 + if self.mode == 9 and self.standalone_op == "dithering" and self.dithering == 0: + raise ValueError("Mode 9 with standalone_op='dithering' requires dithering 1 or 2 (not 0).") + + # Mode 9 image-enhancement methods are standalone transforms; dithering is a separate standalone operation. + if self.mode == 9 and self.standalone_op == "ie_methods" and self.dithering != 0: + raise ValueError("Mode 9 with standalone_op='ie_methods' does not accept dithering. Use standalone_op='dithering' instead.") + + # Image enhancement is a one-shot transform; iterations > 1 has no meaning — silently clamp to 1. + if self.mode == 9 and self.standalone_op == "ie_methods" and self.iterations != 1: + object.__setattr__(self, "iterations", 1) + console_log(msg="Image enhancement is a one-shot transform — iterations forced to 1.", color=Bcolors.WARNING, verbose=self.verbose > 0) # target_hist should match expected images size under as_gray and linear_luminance if self.target_hist is not None and not isinstance(self.target_hist, (str, Path)): diff --git a/src/shinier/SHINIER.py b/src/shinier/SHINIER.py index 1a495ef..125649c 100644 --- a/src/shinier/SHINIER.py +++ b/src/shinier/SHINIER.py @@ -9,7 +9,8 @@ from shinier import ImageDataset, Options, ImageProcessor, REPO_ROOT from shinier.utils import ( Bcolors, console_log, load_np_array, colorize, - print_shinier_header, generate_pydantic_key_value_dict + print_shinier_header, pydantic_model_choices, + IMAGE_ENHANCEMENT_METHODS, ) # Compute repo root as parent of /src/shinier/ @@ -191,6 +192,8 @@ def options_display(opts): if opts.whole_image != 1: types.append('mask') types += ['mode', 'color', 'dithering_memory'] + if opts.mode == 9: + types.append('mode9') if opts.mode == 1: types.append('luminance') if opts.mode in (2, 5, 6, 7, 8): @@ -306,7 +309,7 @@ def SHINIER_CLI(images: Optional[np.ndarray] = None, masks: Optional[np.ndarray] "Histogram + Spectrum", "Spatial frequency + Histogram", "Spectrum + Histogram", - "Dithering only" + "Dithering or HE" ]) opts.mode = mode @@ -316,6 +319,31 @@ def SHINIER_CLI(images: Optional[np.ndarray] = None, masks: Optional[np.ndarray] # --------- Custom Profile --------- if prof == 3: + if mode == 9: + _op = prompt("Standalone operation (mode 9)", default=2, + kind="choice", choices=[ + "Dithering — apply dithering then uint8 cast", + "Image enhancement — transform each image independently", + ]) + if _op == 1: + _dith = prompt("Which dithering method?", default=1, + kind="choice", choices=["Noisy-bit dithering", "Floyd–Steinberg dithering"]) + opts.dithering = _dith + opts.standalone_op = "dithering" + else: + _ie_method_names = list(IMAGE_ENHANCEMENT_METHODS) + _he = prompt("Image enhancement algorithm", default=2, + kind="choice", choices=[ + "Classic HE — global histogram equalization (CDF mapping)", + "TIDHE — tripartite histogram equalization per image", + "RDFHE — recursive dualistic fuzzy histogram equalization", + "NFLDICE — nonlinear fuzzification–linear defuzzification ICE", + "BETCE — bi-entropy curve equalization", + "SFCEF — Sakaguchi-type cost-effective filtering", + ]) + opts.standalone_op = "ie_methods" + opts.ie_methods = _ie_method_names[_he - 1] + as_gray = prompt("Load images as grayscale?", default="No", kind="bool") opts.as_gray = as_gray == 1 linear_luminance = prompt("Are pixel values linearly related to luminance?", default=2, kind='choice', choices=[ @@ -373,17 +401,12 @@ def SHINIER_CLI(images: Optional[np.ndarray] = None, masks: Optional[np.ndarray] opts.rec_standard = rec_standard opts.conserve_memory = prompt("Conserve memory (creates a temporary directory and keep only one image in RAM)?", default='y', kind="bool") - # Dithering dith_choices = ["No dithering", "Noisy-bit dithering", "Floyd–Steinberg dithering"] if mode != 9: dith = prompt("Apply dithering before final uint8 cast?", default=1, kind="choice", choices=dith_choices) opts.dithering = dith - 1 - else: - dith = prompt("Which dithering is going to be applied?", default=1, - kind="choice", choices=dith_choices[1:]) - opts.dithering = dith # Seed now = datetime.now() @@ -413,7 +436,7 @@ def SHINIER_CLI(images: Optional[np.ndarray] = None, masks: Optional[np.ndarray] image_exts = "/".join(f".{ext}" for ext in ACCEPTED_FORMATS) thp1 = prompt("How should the target histogram be defined?", default=1, kind="choice", choices=[ 'Average histogram of input images [default]', - 'Flat histogram a.k.a. `histogram equalization`', + 'Equal target histogram (a.k.a. flat histogram or `histogram equalization`)', 'Derive histogram from an input image file', 'Load histogram from a precomputed array (.npy)' ]) @@ -551,7 +574,12 @@ def main(): fig.savefig(args.save_path, dpi=150) print(f"Processing overview figure saved successfully at: {args.save_path}") else: - plt.show() + import matplotlib + if matplotlib.get_backend().lower() == "agg": + print("Non-interactive backend detected — use --save-path to save the figure to a file.") + plt.close(fig) + else: + plt.show() if __name__ == "__main__": diff --git a/src/shinier/__init__.py b/src/shinier/__init__.py index 6cd1718..ae8c30a 100644 --- a/src/shinier/__init__.py +++ b/src/shinier/__init__.py @@ -55,5 +55,5 @@ from .ImageListIO import ImageListIO from .ImageProcessor import ImageProcessor from .SHINIER import SHINIER_CLI -from .utils import StimulusMasker +from .utils import StimulusMasker, imstats, ImageStats from . import color diff --git a/src/shinier/color/Converter.py b/src/shinier/color/Converter.py index 3e7a65d..bed69ea 100644 --- a/src/shinier/color/Converter.py +++ b/src/shinier/color/Converter.py @@ -29,11 +29,12 @@ RGB_STANDARD = Literal["equal", "rec601", "rec709", "rec2020"] RGB2GRAY_WEIGHTS = { - 'equal': [1/3, 1/3, 1/3], - 'rec601': M_RGB2XYZ_601[1, :], - 'rec709': M_RGB2XYZ_709[1, :], - 'rec2020': M_RGB2XYZ_2020[1, :], + 'equal': np.array([1/3, 1/3, 1/3], dtype=np.float64), + 'rec601': M_RGB2XYZ_601[1, :].astype(np.float64), + 'rec709': M_RGB2XYZ_709[1, :].astype(np.float64), + 'rec2020': M_RGB2XYZ_2020[1, :].astype(np.float64), } +MATLAB_RGB2GRAY_WEIGHTS = np.array([0.298936021293775, 0.587043074451121, 0.114020904255103], dtype=np.float64) for k, v in RGB2GRAY_WEIGHTS.items(): RGB2GRAY_WEIGHTS[k] /= np.sum(v) int2key_mapping = dict(zip(range(1, len(RGB2GRAY_WEIGHTS)+1), RGB2GRAY_WEIGHTS.keys())) @@ -317,7 +318,7 @@ def forward_color_treatment( linear_luminance: bool, as_gray: bool, output_other: Optional[ImageListIO] = None, - conversion_type: Literal['sRGB_to_xyY', 'sRGB_to_lab'] = 'sRGB_to_xyY', + color_conversion: Literal['sRGB_to_xyY', 'sRGB_to_lab'] = 'sRGB_to_xyY', desaturate_chroma_on_low_luminance: bool = False, legacy_mode: bool = False, verbose: bool = False) -> Tuple[ImageListIO, Optional[ImageListIO]]: @@ -350,7 +351,7 @@ def forward_color_treatment( output_other : Optional[ImageListIO] Secondary buffer receiving auxiliary chromatic channels. - conversion_type : Literal["sRGB_to_xyY", "sRGB_to_lab"] + color_conversion : Literal["sRGB_to_xyY", "sRGB_to_lab"] Forward color conversion to apply. desaturate_chroma_on_low_luminance : bool @@ -358,7 +359,12 @@ def forward_color_treatment( chromatic noise from being inflated by luminance manipulation. legacy_mode : bool - If True, uses MATLAB-compatible grayscale conversion behavior. + If True, uses MATLAB-compatible grayscale conversion behavior when + ``as_gray=True`` and ``linear_luminance=False``. This matches the Y + channel of MATLAB's ``rgb2ntsc`` transform. When called from the + SHINIER pipeline, ``legacy_mode=True`` forces both ``as_gray=1`` + and ``linear_luminance=False`` via ``Options``, so the condition is + always satisfied. verbose : bool If True, prints processing messages. @@ -401,7 +407,7 @@ def forward_color_treatment( if as_gray: # Convert to grayscale using simple mean for idx, image in enumerate(input_images): - output_images[idx] = rgb2gray(image, conversion_type="equal", matlab_601=legacy_mode) + output_images[idx] = rgb2gray(image, weighting_standard="equal", matlab_601=legacy_mode) else: for idx, image in enumerate(input_images): if np.issubdtype(image.dtype, np.uint8): @@ -411,7 +417,11 @@ def forward_color_treatment( # --- CASE 2: Color treatment branch ------------------------------------- elif not linear_luminance: for idx, image in enumerate(output_images): - if conversion_type == "sRGB_to_xyY": + if as_gray and legacy_mode: + output_images[idx] = rgb2gray(image, weighting_standard="rec601", matlab_601=True) + continue + + if color_conversion == "sRGB_to_xyY": # Convert from sRGB → xyY (internally handles gamma decoding) srgb_before = image / 255.0 _image = converter.sRGB_to_xyY(srgb_before) @@ -435,7 +445,7 @@ def forward_color_treatment( ) output_other[idx] = xy_after - elif conversion_type == "sRGB_to_lab": + elif color_conversion == "sRGB_to_lab": # Convert from sRGB → Lab (internally handles gamma decoding) _image = converter.sRGB_to_lab(image / 255) @@ -446,7 +456,7 @@ def forward_color_treatment( if as_gray == 0: output_other[idx] = _image[:, :, 1:] else: - raise ValueError(f"Unknown conversion type `{conversion_type}`") + raise ValueError(f"Unknown color conversion `{color_conversion}`") return output_images, output_other else: @@ -460,8 +470,9 @@ def backward_color_treatment( linear_luminance: bool, as_gray: bool, input_other: Optional[ImageListIO] = None, - conversion_type: Literal['xyY_to_sRGB', 'lab_to_sRGB'] = 'xyY_to_sRGB', + color_conversion: Literal['xyY_to_sRGB', 'lab_to_sRGB'] = 'xyY_to_sRGB', gamut_strategy: str = 'clip', + legacy_mode: bool = False, verbose: bool = False) -> ImageListIO: """Apply the backward color-treatment step to an image collection. @@ -489,12 +500,16 @@ def backward_color_treatment( input_other : Optional[ImageListIO] Auxiliary chromatic data required for color reconstruction. - conversion_type : Literal["xyY_to_sRGB", "lab_to_sRGB"] + color_conversion : Literal["xyY_to_sRGB", "lab_to_sRGB"] Backward color conversion to apply. gamut_strategy : str Strategy for repairing out-of-gamut pixels during conversion. + legacy_mode : bool + If True and ``as_gray=True``, keep MATLAB-compatible grayscale + intensities without applying an additional sRGB transfer function. + verbose : bool If True, prints processing messages. @@ -538,7 +553,7 @@ def backward_color_treatment( other = input_other[idx] # Rebuild xyY or Lab: shape (H, W, 3) - if conversion_type == "xyY_to_sRGB": + if color_conversion == "xyY_to_sRGB": # --- Out-of-gamut repair --- Y255, other = gamut_control.apply_image(Y=Y*255, other=other, idx=idx, verbose=verbose) Y = Y255/255.0 @@ -546,15 +561,14 @@ def backward_color_treatment( # Gamma Encode to sRGB output_images[idx] = converter.xyY_to_sRGB(np.dstack((other, Y))) * 255 - elif conversion_type == "lab_to_sRGB": + elif color_conversion == "lab_to_sRGB": # Convert xyY → sRGB (includes linear→gamma) lab = np.dstack([Y*100, other]) output_images[idx] = converter.lab_to_sRGB(lab) * 255 # Output = Grayscale image else: - # Apply gamma encoding (sRGB transfer function) - Yg = converter.linRGB_to_sRGB(im3D(Y)) + Yg = im3D(Y) if legacy_mode else converter.linRGB_to_sRGB(im3D(Y)) # Replicate into 3 channels for display compatibility output_images[idx] = np.dstack([Yg, Yg, Yg]) * 255 @@ -562,21 +576,24 @@ def backward_color_treatment( return output_images -def rgb2gray(image: Union[np.ndarray, Image.Image], conversion_type: RGB_STANDARD = 'equal', matlab_601: bool = False) -> np.ndarray: - """Convert an R'G'B' image to grayscale luma. +def rgb2gray( + image: Union[np.ndarray, Image.Image], + weighting_standard: RGB_STANDARD = 'equal', + matlab_601: bool = False, +) -> np.ndarray: + """Convert an RGB image to grayscale using luma or luminance coefficients. Parameters ---------- image : Union[np.ndarray, Image.Image] - RGB image array with a final channel dimension of 3. The image is - assumed to be gamma-encoded, as with typical sRGB files. + RGB image array with a final channel dimension of 3. - conversion_type : RGB_STANDARD - Luma standard to use: ``"equal"``, ``"rec601"``, ``"rec709"``, or - ``"rec2020"``. + weighting_standard : RGB_STANDARD + Weighting standard: ``"equal"``, ``"rec601"``, ``"rec709"``, + or ``"rec2020"``. matlab_601 : bool - If True and ``conversion_type="rec601"``, uses MATLAB's BT.601 weights. + If True and ``weighting_standard="rec601"``, uses MATLAB's BT.601 weights. Returns ------- @@ -585,23 +602,27 @@ def rgb2gray(image: Union[np.ndarray, Image.Image], conversion_type: RGB_STANDAR Notes ----- - This computes luma (Y') from gamma-encoded components. For physical linear - luminance, first linearize RGB, combine linear-light coefficients, then - re-encode if needed. + Rec. standards use the normalized Y row of SHINIER's RGB-to-XYZ matrices. + These linear-light coefficients are an approximation when applied directly + to gamma-encoded RGB values. + + ``matlab_601=True`` uses MATLAB-compatible ``rgb2gray`` / ``rgb2ntsc`` + Y weights (``0.298936021293775, 0.587043074451121, 0.114020904255103``). + Use :func:`rgb2ntsc_intensity` for an explicit NTSC/YIQ intensity image. """ if isinstance(image, Image.Image): image = np.array(image) elif not isinstance(image, np.ndarray): raise ValueError(f"Invalid image type {type(image)}. Supported values are Image.Image and np.ndarray") - if conversion_type not in RGB2GRAY_WEIGHTS.keys(): - raise ValueError('Conversion type must be either rec709, rec601, rec2020 or equal') + if weighting_standard not in RGB2GRAY_WEIGHTS.keys(): + raise ValueError('Weighting standard must be either rec709, rec601, rec2020 or equal') if image.ndim > 2: - if conversion_type == "equal": + if weighting_standard == "equal": return np.mean(image[..., :3], axis=-1) else: - weights = np.array([0.298936021293775, 0.587043074451121, 0.114020904255103]) if conversion_type == "rec601" and matlab_601 else RGB2GRAY_WEIGHTS[conversion_type] + weights = (MATLAB_RGB2GRAY_WEIGHTS if weighting_standard == "rec601" and matlab_601 else RGB2GRAY_WEIGHTS[weighting_standard]).copy() weights /= np.sum(weights) return image[..., 0] * weights[0] + image[..., 1] * weights[1] + image[..., 2] * weights[2] elif image.ndim == 2: @@ -610,6 +631,30 @@ def rgb2gray(image: Union[np.ndarray, Image.Image], conversion_type: RGB_STANDAR raise ValueError(f"Invalid image dimension {image.shape}. Supported values are >= 2") +def rgb2ntsc_intensity(image: Union[np.ndarray, Image.Image]) -> np.ndarray: + """Convert an RGB image to an NTSC intensity image. + + This is a user-facing helper for workflows that specifically request an + NTSC/YIQ intensity image. It is not used by SHINIER's main color-processing + pipeline. The Y channel of MATLAB's ``rgb2ntsc`` transform is equivalent to + MATLAB-compatible ``rgb2gray`` weights: + ``Y = 0.298936021293775 R + 0.587043074451121 G + 0.114020904255103 B``. + The commonly cited NTSC formula ``0.299 R + 0.587 G + 0.114 B`` is the + rounded form of the same coefficients. + + Parameters + ---------- + image : Union[np.ndarray, Image.Image] + RGB image array or PIL image. + + Returns + ------- + np.ndarray + 2D NTSC intensity image. + """ + return rgb2gray(image, weighting_standard="rec601", matlab_601=True) + + def gray2rgb(image: Union[np.ndarray, Image.Image]) -> np.ndarray: """Convert a grayscale image to RGB. diff --git a/src/shinier/color/GamutControl.py b/src/shinier/color/GamutControl.py index 7531995..5697e28 100644 --- a/src/shinier/color/GamutControl.py +++ b/src/shinier/color/GamutControl.py @@ -331,7 +331,7 @@ def apply_low_Y_desaturation( # Optional chroma-loss metric (lazy import to avoid circular imports) if self.log_low_Y_chroma_loss and verbose: try: - from shinier.color.quantify_chroma_loss import relative_mean_chroma_loss_pct_global_lab + from shinier.color.chroma_eval import mean_chroma_loss_pct_lab xy_before = np.asarray(other, dtype=np.float64) xy_after = np.asarray(other_out, dtype=np.float64) @@ -341,7 +341,7 @@ def apply_low_Y_desaturation( srgb_before = self._converter.xyY_to_sRGB(xyY_before) srgb_after = self._converter.xyY_to_sRGB(xyY_after) - loss_pct, mean_c0, mean_c1 = relative_mean_chroma_loss_pct_global_lab( + loss_pct, mean_c0, mean_c1 = mean_chroma_loss_pct_lab( converter=self._converter, srgb_before_01=srgb_before, srgb_after_01=srgb_after, diff --git a/src/shinier/color/__init__.py b/src/shinier/color/__init__.py index 3805e2a..506ecf3 100644 --- a/src/shinier/color/__init__.py +++ b/src/shinier/color/__init__.py @@ -3,14 +3,16 @@ """ from .Converter import ColorConverter, ColorTreatment, COLOR_STANDARDS, WHITE_D65, M_RGB2XYZ_709, M_RGB2XYZ_2020, \ - M_RGB2XYZ_601, rgb2gray, gray2rgb, RGB_STANDARD, REC_STANDARD, RGB2GRAY_WEIGHTS + M_RGB2XYZ_601, rgb2gray, rgb2ntsc_intensity, gray2rgb, RGB_STANDARD, REC_STANDARD, RGB2GRAY_WEIGHTS, \ + MATLAB_RGB2GRAY_WEIGHTS from .GamutControl import GamutControl __all__ = [ "ColorConverter", "ColorTreatment", - "GamutControl" + "GamutControl", "rgb2gray", + "rgb2ntsc_intensity", "gray2rgb", "COLOR_STANDARDS", "WHITE_D65", @@ -19,5 +21,6 @@ "M_RGB2XYZ_709", "REC_STANDARD", "RGB_STANDARD", - "RGB2GRAY_WEIGHTS" + "RGB2GRAY_WEIGHTS", + "MATLAB_RGB2GRAY_WEIGHTS", ] \ No newline at end of file diff --git a/src/shinier/color/quantify_chroma_loss.py b/src/shinier/color/chroma_eval.py similarity index 88% rename from src/shinier/color/quantify_chroma_loss.py rename to src/shinier/color/chroma_eval.py index 14386d3..42f1b7a 100644 --- a/src/shinier/color/quantify_chroma_loss.py +++ b/src/shinier/color/chroma_eval.py @@ -187,28 +187,41 @@ def _quantize_ab( return a_q, b_q, bins -def relative_mean_chroma_loss_pct_global_lab( +def mean_chroma_loss_pct_lab( *, converter: ColorTreatment, srgb_before_01: np.ndarray, srgb_after_01: np.ndarray, eps: float = 1e-12, ) -> Tuple[float, float, float]: - """Compute global relative mean chroma loss in CIE Lab. + """Compute the global relative mean chroma loss percentage in CIELAB. - The metric is: - 100 * (E[C*_before] - E[C*_after]) / E[C*_before] - where: - C* = sqrt(a*^2 + b*^2) + Both sRGB images are converted to CIE Lab, and chroma is computed for each + 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. Parameters - ---------- converter: ColorTreatment instance (used for sRGB->Lab). - srgb_before_01: Original sRGB image in [0,1], shape (H, W, 3). - srgb_after_01: Processed sRGB image in [0,1], shape (H, W, 3). - eps: Small constant to avoid division by zero. + ---------- + converter : ColorTreatment + Converter used for sRGB-to-Lab conversion. + srgb_before_01 : np.ndarray + Original sRGB image in [0, 1], shape ``(H, W, 3)``. + srgb_after_01 : np.ndarray + Processed sRGB image in [0, 1], shape ``(H, W, 3)``. + eps : float + Small constant to avoid division by zero. Returns - ------- A tuple: (loss_pct, mean_c_before, mean_c_after). + ------- + tuple + ``(loss_pct, mean_c_before, mean_c_after)``. """ lab0 = converter.sRGB_to_lab(srgb_before_01) lab1 = converter.sRGB_to_lab(srgb_after_01) @@ -223,7 +236,7 @@ def relative_mean_chroma_loss_pct_global_lab( return float(loss_pct), mean_c0, mean_c1 -def chroma_information_retention_lab_ab( +def lab_chroma_retention( converter: ColorTreatment, srgb_before_01: np.ndarray, srgb_after_01: np.ndarray, @@ -232,7 +245,54 @@ def chroma_information_retention_lab_ab( nbits_per_axis: int = 8, ab_range: Tuple[float, float] = (-128.0, 127.0), ) -> ChromaInfoRetention: - """Compute chroma information retention (CIR) using sparse MI on quantized Lab (a*,b*).""" + """""Compute Chroma Information Retention (CIR) in quantized CIELAB ab* space. + + Both sRGB images are converted to CIELAB, and their a* and b* components + are jointly quantized into discrete chroma states. The function computes + chroma entropy before and after processing, as well as the mutual + information between the original and processed chroma states. + + Chroma information retention is defined as:: + + CIR = I(X; Y) / H(X) + + where ``X`` represents the original Lab a*b* chroma states, ``Y`` represents + the processed Lab a*b* chroma states, ``I(X; Y)`` is their mutual + information, and ``H(X)`` is the entropy of the original chroma states. + + The retention value is clipped to [0, 1]. Values near 1 indicate that most + of the original chroma information is retained, whereas values near 0 + indicate low retention. + + Parameters + ---------- + converter : ColorTreatment + Converter used to transform sRGB images to CIE Lab. + srgb_before_01 : np.ndarray + Original sRGB image in [0, 1], with shape ``(H, W, 3)``. + srgb_after_01 : np.ndarray + Processed sRGB image in [0, 1], with shape ``(H, W, 3)``. + mask : np.ndarray, optional + Boolean mask with shape ``(H, W)`` restricting the calculation to + selected pixels. + nbits_per_axis : int, optional + Number of quantization bits used independently for the Lab a* and b* + components. + ab_range : tuple of float, optional + Minimum and maximum Lab a* and b* values retained before quantization. + Values outside this range are clipped. + + Returns + ------- + ChromaInfoRetention + Lab chroma entropy before and after processing, mutual information, + and chroma information retention. + + Raises + ------ + ValueError + If ``mask`` does not match the spatial dimensions of the images. + """ lab0 = converter.sRGB_to_lab(srgb_before_01) lab1 = converter.sRGB_to_lab(srgb_after_01) @@ -280,6 +340,7 @@ def chroma_information_retention_lab_ab( # safe: denom > 0 by construction, but guard anyway valid = (pxy > 0.0) & (denom > 0.0) + # MI: mutual information, in bits per pixel mi = float(np.sum(pxy[valid] * np.log2(pxy[valid] / denom[valid]))) cir = float(mi / h_x) @@ -295,13 +356,16 @@ def _validate_srgb_01(srgb: np.ndarray) -> np.ndarray: being clipped to [0..1]). Parameters - ---------- srgb: Image array of shape (H, W, 3). Expected range [0, 1]. + ---------- + srgb: Image array of shape (H, W, 3). Expected range [0, 1]. Returns - ------- Float64 sRGB in [0, 1]. + ------- + Float64 sRGB in [0, 1]. Raises - ------ ValueError: If shape is wrong or values are outside [0,1] by a non-trivial margin. + ------ + ValueError: If shape is wrong or values are outside [0,1] by a non-trivial margin. """ if srgb.ndim != 3 or srgb.shape[-1] != 3: raise ValueError(f"Expected sRGB image shape (H,W,3). Got {srgb.shape}") @@ -457,7 +521,7 @@ def _mean_chroma( return float(np.mean(c)) -def compute_chroma_metrics_for_image( +def chroma_metrics_for_image( srgb_01: np.ndarray, *, y1: float, @@ -642,7 +706,7 @@ def aggregate_chroma_metrics( rows: List[AggregateRow] = [] for y1 in y1_values: metrics = [ - compute_chroma_metrics_for_image( + chroma_metrics_for_image( img, y1=float(y1), rec_standard=rec_standard, @@ -718,7 +782,7 @@ def _save_errorbar_plot( plt.close() -def generate_chroma_loss_report( +def build_chroma_loss_report( images_srgb_01: Sequence[np.ndarray], *, y1_values: Union[np.ndarray, Iterable[float]], @@ -1045,8 +1109,8 @@ def md_table(rows: List[AggregateRow]) -> str: # ============================================================================= # Backward-compatible function (entropy-only) # ============================================================================= - -def chroma_info_loss_bits_per_pixel_vs_y1( + +def chroma_info_loss_bpp( srgb_01: np.ndarray, *, rec_standard: REC_STANDARD = "rec709", @@ -1056,13 +1120,49 @@ def chroma_info_loss_bits_per_pixel_vs_y1( measure_mask: Optional[np.ndarray] = None, low_y_only: bool = False, ) -> List[ChromaInfoLossResult]: - """Compute legacy entropy-only chroma information loss vs ``y1``. + """Compute entropy-based chroma information loss in bits per pixel vs ``y1``. + + For each y1 value, the function applies the current low-luminance + desaturation treatment and measures the resulting chroma information loss. + Chroma information is quantified as the joint Shannon entropy before and + after processing. + + Information loss is defined as:: + + ``loss_bpp = H_before - H_after`` + + A positive value indicates a reduction in chroma information after + processing, whereas a negative value indicates an increase. + + Parameters + ---------- + srgb_01 : np.ndarray + sRGB image in [0, 1], with shape (H, W, 3). + rec_standard : REC_STANDARD + Color standard: "rec601", "rec709", or "rec2020". + y1_values : iterable of float + Low-luminance thresholds to evaluate. + nbits_per_axis : int + Quantization bits per CIE Lab a* and b* axis. + ab_range : tuple of float + Clipping range for CIE Lab a* and b* before quantization. + measure_mask : np.ndarray, optional + Boolean mask with shape (H, W) restricting the measurements. + low_y_only : bool + If True, restrict measurements to pixels where Y_orig < y1. + + Returns + ------- + list of ChromaInfoLossResult + Entropy before and after processing, information loss in bits per pixel, + affected-pixel fraction, and number of measured pixels for each y1 + value. Notes ----- - Input must be float in [0, 1]. This function is strict and raises on uint8-style [0, 255] input. - - For richer metrics, use ``generate_chroma_loss_report``. + - For richer metrics, use ``build_chroma_loss_report``. """ srgb = _validate_srgb_01(srgb_01) diff --git a/src/shinier/utils.py b/src/shinier/utils.py index f622830..da9d94d 100644 --- a/src/shinier/utils.py +++ b/src/shinier/utils.py @@ -270,7 +270,7 @@ def rgb2gray(image): """Replicates MATLAB's rgb2gray function (ITU-R rec601).""" if image.ndim >= 3: from shinier.color.Converter import rgb2gray - return rgb2gray(image=image[..., :3], conversion_type='rec601', matlab_601=True) + return rgb2gray(image=image[..., :3], weighting_standard='rec601', matlab_601=True) else: return image @@ -568,8 +568,12 @@ def _blur(self, image: np.ndarray) -> np.ndarray: return np.apply_along_axis(convolve, 1, image) -def get_field_values_from_pydantic_model(field): - """Return all possible categorical values for a Pydantic field.""" +def pydantic_field_choices(field): + """Return all possible categorical values for a Pydantic model field. + + Handles Literal, Union, bool, and Path annotations. Falls back to the + field default when no categorical type is detected. + """ ann = field.annotation def extract_values(ann_type): @@ -628,13 +632,17 @@ def extract_values(ann_type): return unique_vals -def generate_pydantic_key_value_dict(model_cls): - """Return dict of field → possible values for a Pydantic model.""" +def pydantic_model_choices(model_cls): + """Return (possible_values, default_values) for all fields of a Pydantic model. + + Uses ``pydantic_field_choices`` per field to enumerate categorical options + (Literal, Union, bool). Non-categorical fields fall back to their default. + """ possible_values = {} default_values = {} for name, field in model_cls.model_fields.items(): try: - possible_values[name] = get_field_values_from_pydantic_model(field) + possible_values[name] = pydantic_field_choices(field) default = field.default_factory() if getattr(field, "default_factory", None) is not None else field.default default_values[name] = str(default) if isinstance(default, Path) else default except Exception as e: @@ -653,6 +661,7 @@ def hist_plot( descriptives: bool = False, ax: Optional[plt.Axes] = None, show_normalized_rmse: bool = False, + show_gradient_bar: bool = True, ) -> Tuple[plt.Figure, Tuple[Any, Any]]: """Display a histogram with optional target and descriptive statistics. @@ -693,6 +702,10 @@ def hist_plot( show_normalized_rmse : bool If True, shows normalized RMSE between two normalized histograms. + show_gradient_bar : bool + If True (default), appends a grayscale gradient bar below the histogram + using ``make_axes_locatable``. + Returns ------- Tuple[plt.Figure, Tuple[Any, Any]] @@ -785,19 +798,23 @@ def hist_plot( ax_hist.spines[spine].set_visible(False) # --- grayscale gradient bar directly under the histogram axis --- - divider = make_axes_locatable(ax_hist) - ax_bar = divider.append_axes("bottom", size="4%", pad=0.0) # pad=0.0 to stick to the axis - gradient = np.linspace(0, 1, 256, dtype=np.float64).reshape(1, -1) - ax_bar.imshow(gradient, cmap='gray', aspect='auto', extent=[0, 255, 0, 1]) - ax_bar.set_xlim(ax_hist.get_xlim()) - ax_bar.set_xticks([]) - ax_bar.set_yticks([]) - for spine in ax_bar.spines.values(): - spine.set_visible(False) - xlabel_text = "Pixel intensity" - ax_bar.set_xlabel(xlabel_text, fontname=fontname, labelpad=2) - if ax is None: - fig.subplots_adjust(bottom=0.18) + ax_bar = None + if show_gradient_bar: + divider = make_axes_locatable(ax_hist) + ax_bar = divider.append_axes("bottom", size="4%", pad=0.0) # pad=0.0 to stick to the axis + gradient = np.linspace(0, 1, 256, dtype=np.float64).reshape(1, -1) + ax_bar.imshow(gradient, cmap='gray', aspect='auto', extent=[0, 255, 0, 1]) + ax_bar.set_xlim(ax_hist.get_xlim()) + ax_bar.set_xticks([]) + ax_bar.set_yticks([]) + for spine in ax_bar.spines.values(): + spine.set_visible(False) + xlabel_text = "Pixel intensity" + ax_bar.set_xlabel(xlabel_text, fontname=fontname, labelpad=2) + if ax is None: + fig.subplots_adjust(bottom=0.18) + else: + ax_hist.set_xlabel("Pixel intensity", fontname=fontname) # ----------------- descriptives overlay (μ and ±1σ) ----------------- if descriptives: @@ -966,6 +983,7 @@ def imhist_plot( ax_bar, ax_hist = hist_plot( hist=hist_normalized, bins=bins, + title=title if ax is not None else None, target_hist=target_hist, descriptives=descriptives, ax=ax_hist, @@ -984,7 +1002,7 @@ def imhist_plot( def freq_axis(n: int) -> np.ndarray: - """Compute spatial frequency axis for image spectrum""" + """Compute the spatial-frequency axis for image spectra.""" # MATLAB: # even n: -n/2 : n/2-1 # odd n: -n/2 : n/2-1 (with halves → -2.5,-1.5,...,+1.5 for n=5) @@ -997,7 +1015,7 @@ def freq_axis(n: int) -> np.ndarray: def get_radius_grid(x_size: int, y_size: int, legacy_mode: bool = False) -> np.ndarray: - """Compute the radius grid for rotational average""" + """Compute the radius grid for rotational averages.""" f2 = freq_axis(x_size) # rows f1 = freq_axis(y_size) # cols XX, YY = np.meshgrid(f1, f2) # shape (xs, ys) @@ -1327,7 +1345,7 @@ def im_power_spectrum_plot(im: np.ndarray, with_colorbar: bool = True): if arr.ndim == 3 and arr.shape[2] >= 3: from shinier.color.Converter import rgb2gray # suppose rgb2gray dispo; sinon fais la combinaison manuelle - gray = rgb2gray(arr, conversion_type='rec709').astype(np.float64, copy=False) + gray = rgb2gray(arr, weighting_standard='rec709').astype(np.float64, copy=False) else: gray = arr.astype(np.float64, copy=False) @@ -1651,6 +1669,13 @@ def exact_histogram( ------- tuple[np.ndarray, list] Histogram-specified image and order-accuracy values per channel. + + Notes + ----- + The ``'noise'`` strategy intentionally deviates from SHINE's ``match.m``, which uses + one-sided noise ``rand() * 0.1`` in ``[0, 0.1]``. SHINIER uses centered noise in + ``[-noise_level, +noise_level]`` with an adaptive level, avoiding the systematic + upward bias introduced by MATLAB's one-sided approach. """ # --- Validate and prepare inputs --- L = n_bins if n_bins is not None else None @@ -1730,6 +1755,9 @@ def exact_histogram( console_log(msg=msg, indent_level=1, color=Bcolors.WARNING, verbose=True) noise_level = 0.1 elif tie_strategy == 'noise': + # SHINE's legacy match.m uses rand(size(image)) * 0.1, i.e. one-sided noise in [0, 0.1]. + # We intentionally deviate: centered noise in [-noise_level, +noise_level] with an + # adaptive level is better at breaking ties without introducing a systematic bias. noise_level = tie_breaking_noise_level(image) im_sort = image.astype(np.float64, copy=True) if tie_strategy == 'noise' or hybrid_extra_step: @@ -1994,6 +2022,14 @@ def soft_clip(arr: np.ndarray, ------- np.ndarray Clipped and possibly rescaled array. + + Notes + ----- + Frequency-based matching routinely produces out-of-range values; MATLAB + silently hard-clips them via a ``uint8`` cast. SHINIER intentionally + deviates from this even in ``legacy_mode``: hard clipping distorts the + matched spectrum, so ``soft_clip`` rescales the array to limit clipping to + at most ``max_percent`` of values, preserving the distribution shape. """ def _zero_clip_mean_preserving(arr, a, b): @@ -2902,6 +2938,109 @@ def _coat_check( return msg +def tidhe_hist_plot( + img: np.ndarray, + ax: Optional[plt.Axes] = None, + title: Optional[str] = None, + figsize: tuple = (9, 5), + dpi: int = 100, + show_gradient_bar: bool = True, +) -> Tuple[plt.Figure, plt.Axes]: + """Plot the TIDHE histogram decomposition for a 2D image. + + The histogram is divided at pl_l and pl_u into lower, middle, and + upper sub-histograms containing approximately equal pixel proportions + (P_l, P_m, and P_u). Dashed lines show the clipping level of + each sub-histogram, computed from its mean and median bin counts. The three + intensity ranges are inferred from the partitioning levels. + + + Parameters + ---------- + img : np.ndarray + 2D intensity image. If a 3D array is provided, the first channel is used + as-is; the function does not infer color-space semantics. Values are + rounded and clipped to uint8. + ax : plt.Axes, optional + Existing axes to draw into. If omitted, a new figure is created. + title : str, optional + Plot title. + figsize : tuple + Figure size when ``ax`` is not provided. + dpi : int + Figure DPI when ``ax`` is not provided. + show_gradient_bar : bool + If True (default), appends a grayscale gradient bar. + + Returns + ------- + Tuple[plt.Figure, plt.Axes] + Figure and main histogram axes. + """ + # Convert to uint8 and compute TIDHE parameters + gray_u8 = np.clip(np.round((img[..., 0] if img.ndim == 3 else img).astype(np.float64)), 0, 255).astype(np.uint8) + hist_raw = np.bincount(gray_u8.ravel(), minlength=256).astype(np.float64) + n_pixels = hist_raw.sum() + cdf = np.cumsum(hist_raw) / n_pixels + pl_l = np.intp(np.clip(np.argmin(np.abs(cdf - 1.0 / 3.0)), 0, 253)) + pl_u = np.intp(np.clip(np.argmin(np.abs(cdf - 2.0 / 3.0)), pl_l + 1, 254)) + bands = ((0, pl_l), (pl_l + 1, pl_u), (pl_u + 1, 255)) + names = ("l", "m", "u") + clip_levels_raw = [(np.mean(hist_raw[s:e + 1]) + np.median(hist_raw[s:e + 1])) / 2.0 for s, e in bands] + fractions = [hist_raw[s:e + 1].sum() / n_pixels for s, e in bands] + + # Normalize for hist_plot (probability frequencies) + hist_norm = hist_raw / n_pixels + clip_levels_norm = [cl / n_pixels for cl in clip_levels_raw] + + # Base histogram via hist_plot + fig, (ax_bar, ax_hist) = hist_plot( + hist=hist_norm, + ax=ax, + title=title, + figsize=figsize, + dpi=dpi, + show_gradient_bar=show_gradient_bar, + ) + legend = ax_hist.get_legend() + if legend is not None: + legend.remove() + + # TIDHE-specific overlays + xs = np.arange(256) + y_top = max(hist_norm.max(), max(clip_levels_norm)) * 1.18 + for (s, e), color, cl_norm in zip(bands, ("#d8e7f7", "#ddf0df", "#f8dfcd"), clip_levels_norm): + ax_hist.fill_between(xs[s:e + 1], hist_norm[s:e + 1], color=color, zorder=0) + ax_hist.hlines(cl_norm, s, e, colors="#7e2f8e", lw=1.8, linestyles="--") + ax_hist.vlines([pl_l, pl_u], 0, y_top, colors="#d95319", lw=2.0) + ax_hist.set_ylim(0, y_top) + + # Summary box: p_l/p_u are partition levels; P_l/P_m/P_u are pixel proportions. + summary = "\n".join( + [rf"$p_l={pl_l},\ p_u={pl_u}$"] + + [rf"$P_l={fractions[0]:.3f},\ P_m={fractions[1]:.3f},\ P_u={fractions[2]:.3f}$"] + + [rf"$cl_l={clip_levels_raw[0]:.0f},\ cl_m={clip_levels_raw[1]:.0f},\ cl_u={clip_levels_raw[2]:.0f}$"] + ) + ax_hist.text( + 0.985, 0.955, summary, + transform=ax_hist.transAxes, + ha="right", va="top", fontsize=8, color="0.18", + bbox={"boxstyle": "round,pad=0.28", "facecolor": "white", "edgecolor": "0.82", "alpha": 0.88}, + ) + + # Partitioning level ticks on gradient bar (or ax_hist if no bar) + if ax_bar is not None: + ax_bar.set_xticks([pl_l, pl_u]) + ax_bar.set_xticklabels([rf"$pl_l={pl_l}$", rf"$pl_u={pl_u}$"]) + ax_bar.tick_params(axis="x", pad=2, length=3, labelsize=9) + ax_bar.xaxis.label.set_size(10) + else: + ax_hist.set_xticks([pl_l, pl_u]) + ax_hist.set_xticklabels([rf"$pl_l={pl_l}$", rf"$pl_u={pl_u}$"]) + + return fig, ax_hist + + def show_processing_overview(processor: ImageProcessor, img_idx: int = 0, show_figure: bool = True, show_initial_target: bool = False) -> plt.Figure: """Display before/after images and diagnostics for all processing steps in one figure. @@ -2936,8 +3075,10 @@ def show_processing_overview(processor: ImageProcessor, img_idx: int = 0, show_f >>> fig = show_processing_overview(processor, img_idx=0, show_figure=False) """ - import os, matplotlib - if not show_figure and os.environ.get("DISPLAY", "") == "": + import os, sys, matplotlib + # Use Agg only on headless Linux (no DISPLAY/WAYLAND); macOS and Windows have native GUI backends + _headless = sys.platform.startswith("linux") and not os.environ.get("DISPLAY") and not os.environ.get("WAYLAND_DISPLAY") + if not show_figure and _headless: matplotlib.use("Agg") fontname = 'Arial' @@ -3131,6 +3272,38 @@ def _load_overview_image(path: Optional[Path]) -> np.ndarray: axL.set_title(f"Before – spectrum\n({row_label})", fontsize=10, fontname=fontname, pad=8) axR.set_title(f"After – spectrum\n({row_label})", fontsize=10, fontname=fontname, pad=8) + # ---- Histogram equalization (HE or TIDHE) ---- + elif base_step == "ie_methods": + op = getattr(processor.options, "ie_methods", "tidhe") + if op == "tidhe": + tidhe_hist_plot( + img=processor._initial_buffer[img_idx], + ax=axL, + title="Input histogram and TIDHE partition", + ) + tidhe_hist_plot( + img=processor._final_buffer[img_idx], + ax=axR, + title="Output histogram after TIDHE", + ) + tidhe_ymax = max(axL.get_ylim()[1], axR.get_ylim()[1]) + axL.set_ylim(0, tidhe_ymax) + axR.set_ylim(0, tidhe_ymax) + elif op in ("classic_he", "rdfhe", "betce", "nfldice"): + ce_label = "classic histogram equalization" if op == "classic_he" else op.upper() + _ = imhist_plot( + img=processor._initial_buffer[img_idx], + descriptives=False, + title=f"Before – {ce_label}", + ax=axL, + ) + _ = imhist_plot( + img=processor._final_buffer[img_idx], + descriptives=False, + title=f"After – {ce_label}", + ax=axR, + ) + # ---- Dithering ---- elif base_step == "dithering": # Dithering only affects appearance (already visible in row 1) @@ -3302,8 +3475,8 @@ def ssim_sens(image1: np.ndarray, image2: np.ndarray, data_range: Optional[float all_sens = [] all_mssim = [] for ch in range(C): - X = img1_3D[:, :, ch] - Y = img2_3D[:, :, ch] + X = img1_3D[:, :, ch].astype(np.float64) + Y = img2_3D[:, :, ch].astype(np.float64) # Local means (Gaussian) ux = convolve_2d(X, g1d) @@ -3376,18 +3549,26 @@ def ssim_sens(image1: np.ndarray, image2: np.ndarray, data_range: Optional[float class StepSizeController: """Three-regime adaptive step-size controller for SSIM optimization. - Behavior: - A) If SSIM increases noticeably → accept, keep weight. - B) If SSIM increases only slightly (stall) → restart with larger weight. - C) If SSIM decreases → restart with smaller weight. + Behavior + -------- + A. If SSIM increases noticeably, accept and keep the weight. + B. If SSIM increases only slightly, restart with a larger weight. + C. If SSIM decreases, restart with a smaller weight. Attributes - ---------- gain_up (float): Multiplier when escaping a stall (default=1.3). - gain_down (float): Multiplier when correcting overshoot (default=0.6). - stall_thresh (float): ΔSSIM below which we consider a stall (default=1e-5). - alpha_min (float): Minimum allowed step size. - alpha_max (float): Maximum allowed step size. - max_stall_iter (int): Maximum allowed number of stalled iterations. + ---------- + gain_up : float + Multiplier when escaping a stall. + gain_down : float + Multiplier when correcting overshoot. + stall_thresh : float + SSIM change below which the optimization is considered stalled. + alpha_min : float + Minimum allowed step size. + alpha_max : float + Maximum allowed step size. + max_stall_iter : int + Maximum allowed number of stalled iterations. """ def __init__( @@ -3420,17 +3601,20 @@ def update(self, alpha: float, ssim_new: float, Y_new: np.ndarray, gradient_new: """Update step size and decide whether to restart or accept. Parameters - ---------- alpha (float): Current step-size weight. - ssim_new (float): Current mean SSIM value. - Y_new (np.ndarray): Current image. - gradient_new (np.ndarray): Current gradient from ssim_sens() + ---------- + alpha : float + Current step-size weight. + ssim_new : float + Current mean SSIM value. + Y_new : np.ndarray + Current image. + gradient_new : np.ndarray + Current gradient from ``ssim_sens``. Returns - ------- (alpha_new, Y_next, restart) - alpha_new : Updated step-size weight. - Y_next : Either new image or reverted one. - gradient_next : Either new gradient or reverted one. - restart : Whether a restart is needed. + ------- + tuple + ``(alpha_new, ssim_next, Y_next, gradient_next, restart, done)``. """ restart = False self.restart_reason = " " @@ -3487,6 +3671,70 @@ def update(self, alpha: float, ssim_new: float, Y_new: np.ndarray, gradient_new: return alpha, ssim_new, Y_new, gradient_new, restart, done +def mean_spectrum(magnitudes: Sequence[np.ndarray], shape: tuple) -> np.ndarray: + """Average magnitude spectrum across a list of spectra. + + Parameters + ---------- + magnitudes : Sequence[np.ndarray] + Per-image magnitude spectra to average. + shape : tuple + Shape of a single image, used to allocate the accumulator. A trailing + channel axis is added if missing. + + Returns + ------- + np.ndarray + Mean spectrum as a 3D array ``(H, W, C)``. + """ + acc = im3D(np.zeros(shape)) + for mag in magnitudes: + acc += mag + acc /= len(magnitudes) + return im3D(acc) + + +def corr_rmse_per_item( + items: Sequence[np.ndarray], + target: np.ndarray, + normalize_rmse: bool = False, + mode: Literal["assume [0, 1]", "actual range", "expected range", "histogram"] = "actual range", +) -> Tuple[np.ndarray, np.ndarray]: + """Correlation coefficient and RMSE of each item against a shared target. + + Each item and the target are flattened before comparison. + + Parameters + ---------- + items : Sequence[np.ndarray] + Arrays to compare against ``target`` (e.g. per-image histograms or + spectra). + target : np.ndarray + Reference array compared against every item. + normalize_rmse : bool, optional + If True, use :func:`normalized_rmse` (with ``mode``); otherwise use + :func:`compute_rmse`. + mode : str, optional + Normalization mode forwarded to :func:`normalized_rmse`. + + Returns + ------- + tuple[np.ndarray, np.ndarray] + Correlation coefficients and RMSE values, one per item. + """ + target_flat = target.ravel() + N = len(items) + corr, rms = np.zeros((N,)), np.zeros((N,)) + for idx, item in enumerate(items): + item_flat = item.ravel() + corr[idx] = np.corrcoef(item_flat, target_flat)[0, 1] + if normalize_rmse: + rms[idx] = normalized_rmse(item_flat, target_flat, mode=mode) + else: + rms[idx] = compute_rmse(item_flat, target_flat) + return corr, rms + + def hist_match_validation(images: ImageListIO, binary_masks: List[np.ndarray], target_hist: Optional[np.ndarray] = None, normalize_rmse: bool = False) -> Tuple[np.ndarray, np.ndarray]: """Validate histogram matching with correlation and RMSE metrics. @@ -3526,18 +3774,16 @@ def hist_match_validation(images: ImageListIO, binary_masks: List[np.ndarray], t target_hist = target_hist / (target_hist.sum(axis=0, keepdims=True) + 1e-12) # Compute metric - N = len(initial_hist) - corr, rms = np.zeros((N,)), np.zeros((N,)) - for idx, a_hist in enumerate(initial_hist): - corr[idx] = np.corrcoef(a_hist.ravel(), target_hist.ravel())[0, 1] - if normalize_rmse: - rms[idx] = normalized_rmse(a_hist.ravel(), target_hist.ravel(), mode='histogram') - else: - rms[idx] = compute_rmse(a_hist.ravel(), target_hist.ravel()) - return corr, rms + return corr_rmse_per_item(initial_hist, target_hist, normalize_rmse, mode='histogram') -def sf_match_validation(images: ImageListIO, target_spectrum: Optional[np.ndarray] = None, normalize_rmse: bool = False) -> Tuple[np.ndarray, np.ndarray]: +def sf_match_validation( + images: ImageListIO, + target_spectrum: Optional[np.ndarray] = None, + normalize_rmse: bool = False, + fft_padding_mode: Literal[0, 1, 2, 3] = 0, + fft_padding_value: Union[int, Literal[300]] = 300, +) -> Tuple[np.ndarray, np.ndarray]: """Validate spatial-frequency matching with correlation and RMSE metrics. Validates spectral match between a set of input images by comparing their @@ -3552,6 +3798,12 @@ def sf_match_validation(images: ImageListIO, target_spectrum: Optional[np.ndarra used. normalize_rmse : bool, optional If True, return normalized RMSE. + fft_padding_mode : int, optional + Padding mode for FFT (0 = no padding). Must match the mode used to + produce ``target_spectrum`` so that source and target spectra are at + the same spatial scale. + fft_padding_value : int, optional + Constant padding value when ``fft_padding_mode=3``. Returns ------- @@ -3559,20 +3811,21 @@ def sf_match_validation(images: ImageListIO, target_spectrum: Optional[np.ndarra Correlation coefficients and RMSE values. """ - x_size, y_size = images[0].shape[:2] n_channels = 1 if images[0].ndim == 2 else 3 - # Compute spectra and mean spectrum if required - magnitudes, phases = get_images_spectra(images=images) + # Compute spectra with the same padding as the target + magnitudes, phases = get_images_spectra( + images=images, + fft_padding_mode=fft_padding_mode, + fft_padding_value=fft_padding_value, + ) if target_spectrum is None: - target_spectrum = im3D(np.zeros(images[0].shape)) - for idx, mag in enumerate(magnitudes): - target_spectrum += mag - target_spectrum /= len(magnitudes) + target_spectrum = mean_spectrum(magnitudes, images[0].shape) target_spectrum = im3D(target_spectrum) - # Returns the frequencies of the image, bins range from -0.5f to 0.5f (0.5f is the Nyquist frequency) 1/y_size is the distance between each pixel in the image - r_int = get_radius_grid(x_size=x_size, y_size=y_size) + # Radius grid at the (padded) spectrum scale — must match source magnitude shape + ts_x, ts_y = target_spectrum.shape[:2] + r_int = get_radius_grid(x_size=ts_x, y_size=ts_y) target_rot_avg = [] initial_rot_avg = [] for idx, image in enumerate(images): @@ -3591,18 +3844,16 @@ def sf_match_validation(images: ImageListIO, target_spectrum: Optional[np.ndarra initial_rot_avg.append(np.stack(ira).T) # Compute metrics - N = len(initial_rot_avg) - corr, rms = np.zeros((N,)), np.zeros((N,)) - for idx, ira in enumerate(initial_rot_avg): - corr[idx] = np.corrcoef(ira.ravel(), target_rot_avg.ravel())[0, 1] - if normalize_rmse: - rms[idx] = normalized_rmse(ira.ravel(), target_rot_avg.ravel(), mode='actual range') - else: - rms[idx] = compute_rmse(ira.ravel(), target_rot_avg.ravel()) - return corr, rms + return corr_rmse_per_item(initial_rot_avg, target_rot_avg, normalize_rmse, mode='actual range') -def spec_match_validation(images: ImageListIO, target_spectrum: Optional[np.ndarray] = None, normalize_rmse: bool = False) -> Tuple[np.ndarray, np.ndarray]: +def spec_match_validation( + images: ImageListIO, + target_spectrum: Optional[np.ndarray] = None, + normalize_rmse: bool = False, + fft_padding_mode: Literal[0, 1, 2, 3] = 0, + fft_padding_value: Union[int, Literal[300]] = 300, +) -> Tuple[np.ndarray, np.ndarray]: """Validate Fourier-spectrum matching with correlation and RMSE metrics. Validates spectral matching of input images by comparing the spectra of each @@ -3620,13 +3871,23 @@ def spec_match_validation(images: ImageListIO, target_spectrum: Optional[np.ndar used. normalize_rmse : bool, optional If True, return normalized RMSE. + fft_padding_mode : int, optional + Padding mode for FFT (0 = no padding). Must match the mode used to + produce ``target_spectrum`` so that source and target spectra are at + the same spatial scale. + fft_padding_value : int, optional + Constant padding value when ``fft_padding_mode=3``. Returns ------- tuple[np.ndarray, np.ndarray] Correlation coefficients and RMSE values. """ - magnitudes, phases = get_images_spectra(images=images) + magnitudes, phases = get_images_spectra( + images=images, + fft_padding_mode=fft_padding_mode, + fft_padding_value=fft_padding_value, + ) n_channels = 1 if images[0].ndim == 2 else 3 compute_target_spectrum = target_spectrum is None @@ -3635,24 +3896,10 @@ def spec_match_validation(images: ImageListIO, target_spectrum: Optional[np.ndar if n_channels_ts != n_channels: raise ValueError(f"Target spectrum has {n_channels_ts} channels but images have {n_channels} channels") - target_spectrum = im3D(np.zeros(images[0].shape)) if compute_target_spectrum else im3D(target_spectrum) - if compute_target_spectrum: - for idx, mag in enumerate(magnitudes): - target_spectrum += mag - target_spectrum /= len(magnitudes) - target_spectrum = im3D(target_spectrum) + target_spectrum = mean_spectrum(magnitudes, images[0].shape) if compute_target_spectrum else im3D(target_spectrum) # Compute metric - N = len(magnitudes) - corr, rms = np.zeros((N,)), np.zeros((N,)) - for idx, a_mag in enumerate(magnitudes): - a_mag = im3D(a_mag) - corr[idx] = np.corrcoef(a_mag.ravel(), target_spectrum.ravel())[0, 1] - if normalize_rmse: - rms[idx] = normalized_rmse(a_mag.ravel(), target_spectrum.ravel(), mode='actual range') - else: - rms[idx] = compute_rmse(a_mag.ravel(), target_spectrum.ravel()) - return corr, rms + return corr_rmse_per_item(magnitudes, target_spectrum, normalize_rmse, mode='actual range') def compute_rmse(image1: np.ndarray, image2: np.ndarray, log: bool = False) -> float: @@ -3885,7 +4132,7 @@ def rescale_images255(images: ImageListIO, rescaling_option: Literal[0, 1, 2, 3] rescaling_option : Literal[0, 1, 2, 3], optional Rescaling strategy. ``0`` disables rescaling, ``1`` rescales each image independently, ``2`` uses dataset absolute min/max, and ``3`` uses - dataset average min/max. + dataset average min/max (could necessitate clipping). legacy_mode : bool, optional If True, convert results with MATLAB-compatible uint8 behavior. @@ -3920,7 +4167,8 @@ def rescale_images255(images: ImageListIO, rescaling_option: Literal[0, 1, 2, 3] if rescaling_option == 1: new_image = stretch(new_image) * 255 else: - new_image = (new_image - mn)/(mx - mn) * 255 + # Clipping needed with Option 3. mn/mx are dataset averages, values could exceed [0, 255]. + new_image = np.clip((new_image - mn) / (mx - mn) * 255, 0, 255) images[idx] = MatlabOperators.uint8(new_image) if legacy_mode else new_image return images @@ -4073,6 +4321,405 @@ def imhist(image: np.ndarray, mask: Optional[np.ndarray] = None, n_bins: int = 2 return count +def _validate_uint8_gray(image: np.ndarray, name: str) -> None: + """Validate the uint8 grayscale domain used by image-enhancement methods.""" + if image.ndim != 2: + raise ValueError(f"{name} expects a 2D grayscale image") + if image.dtype != np.uint8: + raise ValueError(f"{name} expects dtype uint8") + + +def _hist256(image: np.ndarray) -> np.ndarray: + """Compute a 256-bin uint8 histogram as float64.""" + return np.bincount(image.ravel(), minlength=256).astype(np.float64) + + +def _apply_uint8_lut(image: np.ndarray, mapping: np.ndarray, legacy_mode: bool) -> np.ndarray: + """Apply a 256-entry lookup table to a uint8 image and cast the result.""" + values = mapping[image] + if legacy_mode: + return MatlabOperators.uint8(values) + return np.clip(np.round(values), 0, 255).astype(np.uint8) + + +def classic_he_gray(image: np.ndarray, legacy_mode: bool = False) -> np.ndarray: + """Apply classic global histogram equalization to a 2D grayscale image. + + Maps each intensity level with a normalized cumulative distribution + function (CDF), so the first non-empty input level maps to 0 and the last + maps to 255. The resulting output histogram is approximately flat over + [0, 255]. + + This is the standard (classic) CDF-based HE — it maximizes contrast but can + over-enhance noise and produce washed-out results on natural images. + For a more controlled alternative see :func:`tidhe_gray`. + + Parameters + ---------- + image : np.ndarray + 2D uint8 grayscale image. + legacy_mode : bool, optional + If True, use MATLAB-compatible rounding and uint8 cast via + :class:`shinier.utils.MatlabOperators`. Default is False. + """ + _validate_uint8_gray(image, "classic_he_gray") + + hist = _hist256(image) + N = np.float64(hist.sum()) + if N == 0: + return image.copy() + + cdf = np.cumsum(hist) + cdf_min = cdf[np.flatnonzero(hist)[0]] + denom = N - cdf_min + if denom == 0: + return image.copy() + + mapping = np.float64(255.0) * (cdf - cdf_min) / denom + return _apply_uint8_lut(image, mapping, legacy_mode) + + +def tidhe_gray(image: np.ndarray, legacy_mode: bool = True) -> np.ndarray: + """Apply TRIPARTITE IMAGE DECOMPOSITION-BASED HISTOGRAM EQUALIZATION (TIDHE). + + TIDHE is a state-of-the-art HE-based contrast-enhancement + algorithm designed for slightly low-contrast and low-contrast grayscale + images. It partitions the image into three intensity sub-images, clips + their histograms to control the enhancement rate, and equalizes each + sub-image independently to limit brightness shifts and artifacts. + + Parameters + ---------- + image : np.ndarray + 2D uint8 grayscale image. + legacy_mode : bool, optional + If True (default), use MATLAB-compatible rounding (round-half-away-from-zero) + for clipping levels and the final uint8 cast, reproducing ``imTIDHE`` exactly. + If False, use standard NumPy rounding and casting. + + References + ---------- + 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 + """ + _validate_uint8_gray(image, "tidhe_gray") + + hist = _hist256(image) + N_pixels = np.float64(hist.sum()) + if N_pixels == 0: + return image.copy() + + cdf = np.cumsum(hist) / N_pixels + mapping = np.arange(256, dtype=np.float64) + + # Eqs. (1)-(2) - pl_l, pl_u : partitioning levels where cdf ≈ 1/3 and ≈ 2/3 + # clipped so each of the three bands has at least one bin + pl_l = np.intp(np.clip(np.argmin(np.abs(cdf - np.float64(1.0 / 3.0))), 0, 253)) + pl_u = np.intp(np.clip(np.argmin(np.abs(cdf - np.float64(2.0 / 3.0))), pl_l + 1, 254)) + + # Produce y_l′, y_m′, y_u′ by equalizing each sub-band — Eqs. (6)-(11) + for start, end in ((0, pl_l), (pl_l + 1, pl_u), (pl_u + 1, 255)): + hist_sub = hist[start:end + 1] + + # IHC clipping level: average of mean and median — Eqs. (3)-(5) + raw_cl = np.float64((hist_sub.mean() + np.median(hist_sub)) / 2.0) + cl = MatlabOperators.round(raw_cl) if legacy_mode else np.round(raw_cl) + + # IHC-clipped sub-histogram + ihc = np.minimum(hist_sub, cl) + sum_ihc = np.float64(ihc.sum()) + + if sum_ihc > 0: + # t_l/t_m/t_u transformation functions — Eqs. (9)-(11) + mapping[start:end + 1] = np.float64(start) + np.float64(end - start) * np.cumsum(ihc) / sum_ihc + + # Eq. (12) - Apply LUT to every pixel and cast to uint8 + return _apply_uint8_lut(image, mapping, legacy_mode) + + +def rdfhe_gray(image: np.ndarray, theta: int = 10, legacy_mode: bool = True) -> np.ndarray: + """Apply RECURSIVE DUALISTIC FUZZY HISTOGRAM EQUALIZATION (RDFHE). + + RDFHE is a state-of-the-art fuzzy HE-based contrast-enhancement + algorithm designed for low-contrast grayscale images. It constructs + a fuzzy image histogram, recursively partitions it into four intensity + sub-histograms, and equalizes each sub-histogram independently to + enhance contrast while limiting brightness shifts and entropy loss. + + Parameters + ---------- + image : np.ndarray + 2D uint8 grayscale image. + theta : int, optional + Fuzziness-related parameter ϑ. Original article uses ``theta=10``. + legacy_mode : bool, optional + If True (default), use MATLAB-compatible uint8 casting for the final + image, reproducing ``imRDFHE`` exactly. + + References + ---------- + 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 + """ + _validate_uint8_gray(image, "rdfhe_gray") + if theta <= 0: + raise ValueError("rdfhe_gray expects theta > 0") + hist = _hist256(image) + + # Equivalent implementation of Eq. (1) - n_tilde(i): fuzzy image histogram (FIH) + # Pixels with the same intensity contribute identically, so Eq. (1) can be + # computed exactly by convolving the ordinary histogram with the triangular + # membership function. + offsets = np.arange(-(theta - 1), theta, dtype=np.float64) + triangular_membership_weights = np.maximum( + 0.0, + 1.0 - np.abs(offsets) / np.float64(theta), + ) + n_tilde = np.convolve(hist, triangular_membership_weights, mode="same") + N_tilde = np.float64(np.sum(n_tilde)) + if N_tilde == 0: + return image.copy() + + # Eqs. (2), (3), (4) - Partition Levels: pl_1 (≈25% FIH mass), pl_2 (≈50%), pl_3 (≈75%) + cdf = np.cumsum(n_tilde) / N_tilde + pl_1 = np.intp(np.clip(np.argmin(np.abs(cdf - 0.25)), 0, 252)) + pl_2 = np.intp(np.clip(np.argmin(np.abs(cdf - 0.50)), pl_1 + 1, 253)) + pl_3 = np.intp(np.clip(np.argmin(np.abs(cdf - 0.75)), pl_2 + 1, 254)) + + mapping = np.arange(256, dtype=np.float64) + for start, end in ((0, pl_1), (pl_1 + 1, pl_2), (pl_2 + 1, pl_3), (pl_3 + 1, 255)): + n_sub = n_tilde[start:end + 1] # <- n_tilde_1 .. n_tilde_4, Eqs (5), (6), (7), (8) + sum_n_sub = np.float64(np.sum(n_sub)) + if sum_n_sub > 0: + # Eq. (9): equalize the sub-histogram onto its own range [start, end]. + # `start` is the branch offset, `end - start` its output width (t_s, t_d in branch 4). + mapping[start:end + 1] = np.float64(start) + np.float64(end - start) * np.cumsum(n_sub) / sum_n_sub + + # x′ - Apply LUT to every pixel and cast to uint8 + return _apply_uint8_lut(image, mapping, legacy_mode) + + +def nfldice_gray( + image: np.ndarray, + b: float = 10.0, + e_l: float = 5.0, + p_l: float = 127.5, + legacy_mode: bool = True, +) -> np.ndarray: + """Apply NONLINEAR FUZZIFICATION–LINEAR DEFUZZIFICATION-BASED IMAGE CONTRAST ENHANCEMENT (NFLDICE). + + NFLDICE is a state-of-the-art fuzzy set-theoretic image-enhancement algorithm + designed for grayscale images. It maps each gray level to a fuzzy membership + value using a nonlinear logistic fuzzifier, then converts the membership + value back to the intensity range using a linear defuzzifier. Because the + transformation is monotonic and applied independently to each pixel, it + enhances contrast while preserving spatial structure and fine details. + + Parameters + ---------- + image : np.ndarray + 2D uint8 grayscale image. + b : float, optional + Base ``B`` of the nonlinear fuzzifier. Controls the enhancement level and + should be larger than two. Article uses ``b=10`` (default). + e_l : float, optional + Exponent gain ``E_l`` of the nonlinear fuzzifier. Controls the enhancement + level and should be larger than two. Article uses ``e_l=5`` (default). + p_l : float, optional + Gray level ``P_l`` assigned a membership value of 0.5 (the fuzzifier's + inflection point). Article uses ``p_l=127.5`` (default). + legacy_mode : bool, optional + If True (default), use MATLAB-compatible uint8 casting for the final + image, reproducing ``imNFLDICE`` exactly. If False, use standard NumPy + rounding and casting. + + References + ---------- + Rahman, H. (2025). A time-efficient and effective image contrast enhancement + technique using fuzzification and defuzzification. In M. Mahmud, M. S. Kaiser, + A. Bandyopadhyay, K. Ray, & S. A. Mamun (Eds.), Trends in electronics and health + informatics (Lecture Notes in Networks and Systems, Vol. 1034, pp. 45–58). + Springer Nature Singapore. https://doi.org/10.1007/978-981-97-3937-0_4 + """ + _validate_uint8_gray(image, "nfldice_gray") + if b <= 1.0: + raise ValueError("nfldice_gray expects b > 1") + if e_l <= 0.0: + raise ValueError("nfldice_gray expects e_l > 0") + + L = 256 + levels = np.arange(L, dtype=np.float64) + + # Eqs. (1)-(2) - nonlinear fuzzifier F_X: logistic membership in (0, 1). + # P_l is the gray level with membership 0.5; B and E_l set the steepness. + membership = 1.0 / (1.0 + b ** (-e_l * ((levels - p_l) / np.float64(L - 1)))) + + # Eqs. (3)-(4) - linear defuzzifier D_L: rescale membership to [0, L-1]. + mapping = membership * np.float64(L - 1) + + # Apply LUT to every pixel and cast to uint8 + return _apply_uint8_lut(image, mapping, legacy_mode) + + +def betce_gray(image: np.ndarray, legacy_mode: bool = True) -> np.ndarray: + """Apply BI-ENTROPY CURVE EQUALIZATION (BETCE). + + BETCE is a state-of-the-art curve-based image-enhancement 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 to enhance contrast while limiting + mean-brightness shifts and entropy loss. + + Parameters + ---------- + image : np.ndarray + 2D uint8 grayscale image. + legacy_mode : bool, optional + If True (default), use MATLAB-compatible rounding and uint8 casting for + the final image, reproducing ``imBETCE`` exactly. + + References + ---------- + 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 + """ + _validate_uint8_gray(image, "betce_gray") + + L = 256 + hist = _hist256(image) + N_pixels = np.float64(hist.sum()) + if N_pixels == 0: + return image.copy() + prob = hist / N_pixels + + # Eq. (1) - ET(i): entropy curve (ETC) + ET = np.zeros(L, dtype=np.float64) + nz = prob > 0 + ET[nz] = -prob[nz] * np.log2(prob[nz]) + sum_ET = np.float64(np.sum(ET)) + if sum_ET == 0: + return image.copy() + + # Eq. (2) - pl: ETC weighted arithmetic mean + levels = np.arange(L, dtype=np.float64) + raw_pl = np.sum(levels * ET) / sum_ET + pl = np.intp(MatlabOperators.round(raw_pl) if legacy_mode else np.round(raw_pl)) + pl = np.intp(np.clip(pl, 0, L - 2)) + + mapping = levels.copy() + for start, end in ((0, pl), (pl + 1, L - 1)): + ET_sub = ET[start:end + 1] # ET_l and ET_u, Eqs. (3), (4) + sum_ET_sub = np.float64(np.sum(ET_sub)) + if sum_ET_sub > 0: + # Eq. (5) - equalize each sub-ETC onto its own brightness range + mapping[start:end + 1] = np.float64(start) + np.float64(end - start) * np.cumsum(ET_sub) / sum_ET_sub + + return _apply_uint8_lut(image, mapping, legacy_mode) + + +def sfcef_gray(image: np.ndarray, t: float = 0.5, legacy_mode: bool = True) -> np.ndarray: + """Apply SAKAGUCHI-TYPE FUNCTION-BASED COST-EFFECTIVE FILTERING (SFCEF). + + SFCEF is a state-of-the-art filtering-based image-enhancement algorithm + designed for low-light grayscale images. It computes two filter coefficients + from coefficient bounds of a Sakaguchi/Gegenbauer geometric function class, + then convolves the image with a single 3x3 filter. + + Parameters + ---------- + image : np.ndarray + 2D uint8 grayscale image. + t : float, optional + Fixed parameter of the geometric function class. The article uses + ``t=0.5`` for low-light images (default) and ``t=-1`` for low-contrast + images. + legacy_mode : bool, optional + If True (default), use MATLAB-compatible uint8 casting. Exact MATLAB + parity is still not guaranteed because MATLAB ``filter2`` may differ + from NumPy by one gray level around rounding ties. + + Notes + ----- + Python ``sfcef_gray(..., legacy_mode=True)`` was compared against MATLAB + ``imSFCEF`` on the 500 low-light images of the LOL dataset + (Kaggle mirror: https://www.kaggle.com/datasets/soumikrakshit/lol-dataset). + MATLAB outputs were generated as ``imread -> rgb2gray -> imSFCEF`` and Python + outputs used MATLAB-compatible grayscale conversion. The goal of this + comparison was to quantify numerical agreement, not to establish that one + version is better than the other. No pixel differed by more than one gray + level, and 0.45075667% of pixels differed by exactly one gray level. + + Mean metric values (bold = better; max% = worst-case per-image gap vs MATLAB): + + - AMBE (↓): Py **42.7434** / MATLAB 42.7444 [max 0.0080%] + - MSSIM (↑): Py **0.6633** / MATLAB 0.6633 [max 0.0191%] + - PSNR (↑): Py **14.7413** / MATLAB 14.7410 [max 0.0145%] + - BP2BPSIM (↑): Py 0.5372 / MATLAB **0.5372** [max 0.0242%] + - CI (↑): Py 58.2189 / MATLAB **58.2194** [max 0.0300%] + - Entropy (↑): Py 6.8514 / MATLAB **6.8528** [max 0.0638%] + + All differences are small enough to be treated as numerical implementation + differences. + + References + ---------- + 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 + """ + _validate_uint8_gray(image, "sfcef_gray") + + phi = np.float64(0.5) + x = np.float64(1.0) + u2 = np.float64(1.0 + t) + u3 = np.float64(1.0 + t + t * t) + + # Eqs. (1), (2), (3) - coefficient bounds a1, a2, a3 of G_S(phi) + a1 = np.float64(1.0) + a2 = (2.0 * phi * x) ** np.float64(3.0 / 2.0) / np.sqrt( + phi * (2.0 - u2) ** 2.0 + - 2.0 * x ** 2.0 * ( + (phi + phi ** 2.0) * (2.0 - u2) ** 2.0 + - 2.0 * phi ** 2.0 * ((3.0 - u3) - (2.0 * u2 - u2 ** 2.0)) + ) + ) + a3 = (2.0 * phi * x / (2.0 - u2)) ** 2.0 + (2.0 * phi * x) / (3.0 - u3) + + # Eqs. (4)-(5) with d1=d3=1/8, d2=d4=d6=0, d5=1 + c1 = (a1 + a3) / 8.0 + c2 = a2 + + # Steps 5-6 - 3x3 filter: c1 around the center and c2 at center + padded = np.pad(image.astype(np.float64), 1, mode="constant") + neighbours = ( + padded[:-2, :-2] + padded[:-2, 1:-1] + padded[:-2, 2:] + + padded[1:-1, :-2] + padded[1:-1, 2:] + + padded[2:, :-2] + padded[2:, 1:-1] + padded[2:, 2:] + ) + filtered = c1 * neighbours + c2 * padded[1:-1, 1:-1] + + if legacy_mode: + return MatlabOperators.uint8(filtered) + return np.clip(np.round(filtered), 0, 255).astype(np.uint8) + + +IMAGE_ENHANCEMENT_METHODS = { + "classic_he": {"fn": classic_he_gray, "label": "Classic HE"}, + "tidhe": {"fn": tidhe_gray, "label": "TIDHE"}, + "rdfhe": {"fn": rdfhe_gray, "label": "RDFHE"}, + "nfldice": {"fn": nfldice_gray, "label": "NFLDICE"}, + "betce": {"fn": betce_gray, "label": "BETCE"}, + "sfcef": {"fn": sfcef_gray, "label": "SFCEF"}, +} + + def rounded_target_hist(target_hist: np.ndarray, n_pixels: int) -> np.ndarray: """Convert an ideal target histogram into realizable probabilities. @@ -4219,3 +4866,161 @@ def normalize_hist(a_hist): return average, hist_list else: return average + + +@dataclass +class ImageStats: + """Per-image and dataset-level luminance statistics. Mirrors the output of SHINE ``imstats``.""" + + #: ``(n_bins, n_images)`` array of histograms for each image. + hist_mat: np.ndarray + #: ``(n_images,)`` array of mean luminance for each image. + mean_vec: np.ndarray + #: ``(n_images,)`` array of luminance standard deviation for each image. + std_vec: np.ndarray + #: ``(n_bins,)`` array of the average histogram across all images. + mean_hist: np.ndarray + #: Scalar mean luminance across all images. + mean_lum: float + #: Scalar mean standard deviation across all images. + mean_std: float + + +def imstats( + images: Union["ImageListIO", List[np.ndarray]], + mask: Optional[Union[np.ndarray, List[np.ndarray]]] = None, + n_bins: int = 256, +) -> "ImageStats": + """Compute per-image and dataset-level luminance statistics (mirrors SHINE ``imstats``).""" + imgs = list(images) + n = len(imgs) + if n == 0: + raise ValueError("images must not be empty.") + + masks = [None] * n if mask is None else [mask] * n if isinstance(mask, np.ndarray) else list(mask) + if len(masks) != n: + raise ValueError("len(mask) must equal the number of images.") + + mean_vec = np.zeros(n) + std_vec = np.zeros(n) + hist_mat = np.zeros((n_bins, n)) + + for i, im in enumerate(imgs): + im1 = MatlabOperators.rgb2gray(im) if im.ndim == 3 else im + bm = np.ones(im1.shape, dtype=bool) if masks[i] is None else np.asarray(masks[i], dtype=bool) + px = im1[bm] + mean_vec[i] = MatlabOperators.mean2(px) + std_vec[i] = MatlabOperators.std2(px) + hist_mat[:, i] = imhist(px.reshape(-1, 1), n_bins=n_bins)[:, 0] + + return ImageStats( + hist_mat=hist_mat, + mean_vec=mean_vec, + std_vec=std_vec, + mean_hist=hist_mat.mean(axis=1), + mean_lum=float(mean_vec.mean()), + mean_std=float(std_vec.mean()), + ) + + +def _check_same_shape(image1: np.ndarray, image2: np.ndarray, metric_name: str) -> None: + if image1.shape != image2.shape: + raise ValueError(f"{metric_name} expects images with the same shape.") + + +def compute_ambe(reference: np.ndarray, enhanced: np.ndarray) -> float: + """Absolute mean brightness error between two images. + + Lower values indicate better mean-brightness preservation. + """ + _check_same_shape(reference, enhanced, "compute_ambe") + ref = im3D(reference).astype(np.float64) + enh = im3D(enhanced).astype(np.float64) + return np.float64(np.mean(np.abs(np.mean(ref, axis=(0, 1)) - np.mean(enh, axis=(0, 1))))) + + +def compute_contrast_improvement(image: np.ndarray, n_bins: int = 256) -> float: + """Contrast Improvement (CI) metric. + + Computes the histogram-weighted standard deviation of pixel intensities. + Higher values indicate a more spread-out intensity distribution. + """ + arr = im3D(image).astype(np.float64) + ci_vals = [] + for c in range(arr.shape[2]): + hist, _ = np.histogram(arr[:, :, c], bins=n_bins, range=(0, n_bins)) + total = hist.sum() + if total == 0: + ci_vals.append(0.0) + continue + px = hist / total + bins = np.arange(n_bins, dtype=np.float64) + mu = np.float64(np.sum(bins * px)) + ci = np.float64(np.sqrt(np.sum((bins - mu) ** 2 * px))) + ci_vals.append(ci) + return np.float64(np.mean(ci_vals)) + + +def compute_image_entropy(image: np.ndarray, n_bins: int = 256) -> float: + """Shannon entropy of the image intensity distribution. + + Higher values indicate greater intensity diversity. + """ + arr = im3D(image).astype(np.float64) + ent_vals = [] + for c in range(arr.shape[2]): + hist, _ = np.histogram(arr[:, :, c], bins=n_bins, range=(0, n_bins)) + total = hist.sum() + if total == 0: + ent_vals.append(0.0) + continue + px = hist / total + nonzero = px > 0 + ent = np.float64(-np.sum(px[nonzero] * np.log2(px[nonzero]))) + ent_vals.append(ent) + return np.float64(np.mean(ent_vals)) + + +def compute_mssim(reference: np.ndarray, enhanced: np.ndarray, data_range: Optional[float] = None) -> float: + """Mean Structural Similarity Index (MSSIM) between two images. + + Wrapper around :func:`ssim_sens` returning only the per-channel mean SSIM + averaged over all channels. Higher is better (maximum is 1). + """ + _check_same_shape(reference, enhanced, "compute_mssim") + _, ssim_vals = ssim_sens(reference, enhanced, data_range=data_range) + return np.float64(np.mean(ssim_vals)) + + +def compute_psnr(reference: np.ndarray, enhanced: np.ndarray, data_range: float = 255) -> float: + """Peak signal-to-noise ratio between two images. + + Higher values indicate lower reconstruction error. Identical images return + ``np.inf``. + """ + _check_same_shape(reference, enhanced, "compute_psnr") + ref = im3D(reference).astype(np.float64) + enh = im3D(enhanced).astype(np.float64) + values = [] + for ch in range(ref.shape[2]): + mse = np.mean((ref[:, :, ch] - enh[:, :, ch]) ** 2) + values.append(np.inf if mse == 0 else 20.0 * np.log10(np.float64(data_range)) - 10.0 * np.log10(mse)) + return np.float64(np.mean(values)) + + +def compute_bp2bpsim(reference: np.ndarray, enhanced: np.ndarray, n_bits: int = 8) -> float: + """Bit-plane to bit-plane similarity between two images. + + Values range from 0 to 1. Higher values indicate more matching bits across + corresponding pixels and channels. + """ + _check_same_shape(reference, enhanced, "compute_bp2bpsim") + if n_bits < 1 or n_bits > 8: + raise ValueError("compute_bp2bpsim expects n_bits between 1 and 8.") + + ref = np.clip(np.round(reference), 0, 255).astype(np.uint8) + enh = np.clip(np.round(enhanced), 0, 255).astype(np.uint8) + bit_matches = 0 + for bit in range(n_bits): + bit_matches += np.count_nonzero(((ref >> bit) & 1) == ((enh >> bit) & 1)) + return np.float64(bit_matches / (ref.size * n_bits)) diff --git a/tests/README.md b/tests/README.md index 5598099..48ff02b 100644 --- a/tests/README.md +++ b/tests/README.md @@ -58,30 +58,44 @@ pytest -n 4 -s -m unit_tests --- -## ⚡ Validation Tests with Shards +## ⚡ Validation Tests — Coverage Modes -Shards divide the exhaustive test space across multiple processes. +`ImageProcessor_validation_test.py` supports three coverage modes selected via `COVERAGE_MODE`: + +| Mode | Description | Typical shards | +|---------------|---------------------------------------------------------------|----------------| +| `sampled` | Random sample of the full parameter space (default) | 1–10 | +| `pruned` | Full Cartesian product over a reduced parameter set | 10–20 | +| `exhaustive` | Full Cartesian product over every possible combination | 100+ | + +**Pruned mode** collapses parameters whose values share the same code path (e.g. all `ie_methods` go through the same dispatch, `rec_standard` values differ only in RGB→gray weights). +Every collapsed value is independently covered by `tests/unit_tests/PrunedCoverage_test.py`, which verifies that each excluded value runs without error and produces finite output — making pruned-mode results safely extrapolatable. ### Environment variables -| Variable | Description | Default | -| ------------------ | ------------------------------------------ | ------- | -| `SHARDS` | Total number of shards | `1` | -| `SHARD_INDEX` | Index of current shard (0-based) | `0` | -| `SHOW_PROGRESS` | Enable tqdm progress bars | `0` | -| `PERCENT_SAMPLED` | Fraction of combinations to sample | `1` | -| `DUMP_FILE_FORMAT` | Format for failure dumps (`json` or `pkl`) | `json` | -| `START_AT` | Resume testing from given combo index | `0` | +| Variable | Description | Default | +| ------------------ | -------------------------------------------------------- | --------- | +| `COVERAGE_MODE` | `sampled` / `pruned` / `exhaustive` | `sampled` | +| `SHARDS` | Total number of shards | `1` | +| `SHARD_INDEX` | Index of current shard (0-based) | `0` | +| `SHOW_PROGRESS` | Enable tqdm progress bars | `0` | +| `PERCENT_SAMPLED` | Fraction of combinations per shard (sampled mode only) | `0.001` | +| `DUMP_FILE_FORMAT` | Format for failure dumps (`json` or `pkl`) | `json` | +| `START_AT` | Resume exhaustive/pruned from given combo index | `0` | +| `RESTART` | Re-initialize the hash registry before starting | `false` | + +Hash keys are namespaced by `COVERAGE_MODE` (e.g. `exhaustive:5` vs `pruned:5`) so separate runs never collide in the shared `tests/hash_registry.db`. --- -## 🛠️ Running Shards in Parallel +## 🛠️ Running Shards Locally (GNU parallel) Use **GNU parallel** to distribute shards across CPU cores: ```bash parallel --ungroup --jobs 8 \ - 'PYTHONUNBUFFERED=1 SHOW_PROGRESS=1 DUMP_FILE_FORMAT=pkl SHARDS=8 SHARD_INDEX={} pytest -s -m validation_tests' ::: 0 1 2 3 4 5 6 7 + 'COVERAGE_MODE=sampled SHARDS=8 SHARD_INDEX={} SHOW_PROGRESS=1 DUMP_FILE_FORMAT=pkl \ + pytest -s -m validation_tests' ::: 0 1 2 3 4 5 6 7 ``` > 🔹 `--ungroup` allows live tqdm updates in real time.\ @@ -89,6 +103,49 @@ parallel --ungroup --jobs 8 \ --- +## 🖥️ 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`, +`COVERAGE_MODE`, `PERCENT_SAMPLED`, and `--time` to your needs: + +```bash +#!/bin/bash +#SBATCH --job-name=shinier_validation +#SBATCH --array=0- # e.g. 0-99 for 100 shards +#SBATCH --time=48:00:00 # 48 h for exhaustive; 12 h for pruned; 2 h for sampled +#SBATCH --mem=8G +#SBATCH --cpus-per-task=1 +#SBATCH --output=logs/val_%A_%a.out +#SBATCH --error=logs/val_%A_%a.err + +set -euo pipefail +cd "$SLURM_SUBMIT_DIR" +source ....../venv/bin/activate + +COVERAGE_MODE=exhaustive \ # exhaustive | pruned | sampled +SHARDS=$SLURM_ARRAY_TASK_COUNT \ +SHARD_INDEX=$SLURM_ARRAY_TASK_ID \ +START_AT=0 \ # set > 0 to resume after a crash +SHOW_PROGRESS=0 \ +DUMP_FILE_FORMAT=pkl \ +python -m pytest tests/validation_tests/ImageProcessor_validation_test.py \ + -m validation_tests -v --tb=short +``` + +Submit with: +```bash +mkdir -p logs +sbatch run_validation.sh +``` + +To resume exhaustive mode from combo `1180676` after a crash: +```bash +START_AT=1180676 sbatch run_validation.sh +``` + +--- + ## 🔎 Debugging Tests ### Drop into debugger on failure @@ -111,6 +168,28 @@ pytest -m validation_tests -vv -s --tb=long --- +## 🟡 Hard vs Soft Failures in Validation Tests + +`ImageProcessor_validation_test.py` distinguishes two levels of failure: + +| Level | When it fires | Test outcome | +|-------|--------------|--------------| +| **Hard failure** | Unexpected exception or internal validation flag | Test **fails** — something is definitely broken | +| **Soft failure** | Optimization anomaly that is expected in some combos | Test **passes** — issue is printed and dumped, but doesn't block the run | + +### Hard failures (call `_dump_and_fail`) +- Any `proc.validation` entry with `valid_result=False` (e.g. histogram RMSE didn't converge to 0 in single-objective modes). +- RMSE more than doubled or increased > 0.1 in absolute terms in composite modes (5–8). +- SSIM `final` check: the actual output Y has lower SSIM than the first optimization iteration — the optimizer went backward overall. + +### Soft failures (collected in `combo_warnings`, dumped but don't fail) +- RMSE regressed strictly but by a small amount — expected in composite modes (spec_match undoes histogram progress in mode 6, etc.) or when `hist_optim=True` explicitly trades RMSE for SSIM. +- SSIM `sub_iter` check: the last *pre-rollback proposal* of the SSIM optimizer had lower SSIM than the first — a known artifact of the rollback mechanism (the final output may still be fine). + +Soft failures produce a `.pkl` / `.json` dump in `tests/assets/tmp/` for inspection, and are printed as a summary at the end of the test run. + +--- + ## 🔄 Resume From a Given Combo If a bug occurs at combo 21,600 (from tqdm output): @@ -122,6 +201,70 @@ pytest -m validation_tests -vv -s --maxfail=1 --pdb --- +## Image Enhancement MATLAB Reference Hashes + +Each image-enhancement algorithm in `ie_methods` is validated pixel-strictly against its free +MATLAB reference implementation. The test file is +`tests/validation_tests/ImageEnhancement_validation_test.py`. + +**Approach:** reference outputs were generated once from +`tests/assets/SAMPLE_512X512/*.png`, converted to raw `uint8` bytes, and stored +as SHA256 hashes in `tests/assets/image_enhancement_matlab_sha256.json`. +The JSON stores the shared `dtype`, `shape`, and ordered `images` list once; +each algorithm stores only the SHA256 list in that same image order. The Python +test recomputes each output and compares shape, dtype, and hash. + +Full citations are in the *Implemented algorithms* section of +`documentation/documentation.md`. + +| Algorithm | MATLAB fn | Parameters | JSON key | MATLAB code DOI | +|-----------|-------------|-----------------------------|-----------|------------------------------------| +| TIDHE | `imTIDHE` | defaults | `tidhe` | 10.13140/RG.2.2.22946.70088 | +| RDFHE | `imRDFHE` | `p=10` | `rdfhe` | 10.13140/RG.2.2.23921.34408 | +| NFLDICE | `imNFLDICE` | `B=10, E_l=5, P_l=127.5` | `nfldice` | 10.13140/RG.2.2.14716.51849 | +| BETCE | `imBETCE` | defaults | `betce` | 10.13140/RG.2.2.14319.09126 | +| SFCEF | `imSFCEF` | `t=0.5` | `sfcef` | 10.13140/RG.2.2.16448.44807 | + +TIDHE, RDFHE, NFLDICE, and BETCE use exact SHA256 validation against MATLAB +reference hashes. + +SFCEF uses a different validation strategy because MATLAB's `imSFCEF` relies on +`filter2`. MATLAB delegates this convolution to platform libraries that can use +Fused Multiply-Add (FMA) instructions, producing float64 values that differ from +NumPy by a few units in the last place before final `uint8` rounding (see +*MATLAB vs Python Differences §5* in `documentation/documentation.md`). Because +of this, exact SHA256 equality is not a stable requirement for SFCEF. + +The strict SFCEF validation therefore loads MATLAB `imSFCEF` reference images +from `tests/assets/sfcef_matlab_reference/` and asserts that Python +`sfcef_gray(..., legacy_mode=True)` differs by **at most 1 gray level per +pixel** (`max_diff ≤ 1`). The JSON still stores SFCEF MATLAB SHA256 hashes as a +reference record, but SFCEF's pass/fail test uses the pixel-difference bound +instead of exact hash equality. + +We also ran a broader MATLAB-vs-Python comparison on the 500 low-light images of +the LOL dataset +(`https://www.kaggle.com/datasets/soumikrakshit/lol-dataset`). MATLAB outputs +were generated as `imread -> rgb2gray -> imSFCEF`; Python outputs used +MATLAB-compatible grayscale conversion followed by `sfcef_gray(..., +legacy_mode=True)`. No pixel differed by more than one gray level, and +0.45075667% of pixels differed by exactly one gray level. The goal of this +larger comparison was not to claim that one implementation is better, but to +verify that the numerical differences are negligible and do not systematically +favor either implementation. + +Mean metric values on LOL were nearly identical: AMBE against the paired LOL +high image was 42.74335064 for Python and 42.74440659 for MATLAB; MSSIM was +0.66327749 for Python and 0.66325624 for MATLAB; PSNR was 14.74127390 for +Python and 14.74096545 for MATLAB; BP2BPSIM was 0.53723447 for Python and +0.53723595 for MATLAB; CI was 58.21891968 for Python and 58.21938888 for +MATLAB; entropy was 6.85142815 for Python and 6.85281768 for MATLAB. Maximum +absolute metric differences, expressed as percentages of the MATLAB values, +were AMBE 0.007961%, MSSIM 0.019079%, PSNR 0.014468%, BP2BPSIM 0.024168%, CI +0.030034%, and entropy 0.063809%. + +--- + ## 🤍 Replay a Dumped Failure To reproduce a failed validation test: @@ -152,10 +295,12 @@ This will rebuild the same `Options`, reload selected images, and re-run the fai ## 🔧 Example Workflow -1. Run exhaustive validation tests across 8 shards: +1. Run sampled validation tests locally across 8 cores: ```bash - parallel --ungroup --jobs 8 'SHOW_PROGRESS=1 DUMP_FILE_FORMAT=pkl SHARDS=8 SHARD_INDEX={} pytest -m validation_tests -x -s' ::: 0 1 2 3 4 5 6 7 + parallel --ungroup --jobs 8 \ + 'COVERAGE_MODE=sampled SHOW_PROGRESS=1 DUMP_FILE_FORMAT=pkl SHARDS=8 SHARD_INDEX={} \ + pytest -m validation_tests -s' ::: 0 1 2 3 4 5 6 7 ``` 2. Inspect failures: diff --git a/tests/assets/image_enhancement_matlab_sha256.json b/tests/assets/image_enhancement_matlab_sha256.json new file mode 100644 index 0000000..24186e9 --- /dev/null +++ b/tests/assets/image_enhancement_matlab_sha256.json @@ -0,0 +1,86 @@ +{ + "dtype": "uint8", + "shape": [ + 512, + 512 + ], + "images": [ + "test_face1.png", + "test_face2.png", + "test_face3.png", + "test_face4.png", + "test_face5.png", + "test_scene1.png", + "test_scene2.png", + "test_scene3.png", + "test_scene4.png", + "test_scene5.png" + ], + "sha256": { + "tidhe": [ + "8bc0f1c3ab413a42310a3a7f7f22bd5c9d029182205b5fc69abd7538c1d5b09e", + "e59715d62c6f849d41ed0e060ede84df21b923a25ad1d7f3eece61d337c302cc", + "09489df99428b568778a6f98a66d65636b51d857b4a3f19699e9ecf7c7749717", + "043429ee09df6aafc5c76d9ca02751f5fb20169d5fe48289d0c8e69fec53895f", + "fb8dd8fa84571f7bf8ca7c0e17536d41c57c78c8a7820dd6ffc16b7be0803b78", + "a8912207b70e7349d826e0f42875a7beafbd74e5635c18b69ea72f856ca8a51a", + "b117c0d6fa7a54b3f0eb3a3ac637ddaaa67b2db975ce4413a19423c4c33f4dbb", + "8c10e54b3d893e4f1f641c0da8a4fd108a22a529a1f211350a3e51ca31cc70a9", + "02689dc3b3b52af64824704b88a6ca9f98ebace541b51800215502d04beb8772", + "64599bbf3b4df2d41c3a98f75d34ca0a984673266744c9ffa0ae91b4d34aba36" + ], + "rdfhe": [ + "7ed2f181c3058b59666f167a9e8d70bca395dd2e06e23eef7af2be1635359860", + "5615a9d9dd4a1ffa7965ee6cb58e0b4e1b20b97f86a1b55ca3afb6eca1b25705", + "1a908f53eb274c4f11c2da1d3031e0503ea343544bbec12bd9fc4557524dc182", + "9d64110f1f7ea0d30afaaffcd58da7b9ef15de0be3902f87bd5bab5339e76947", + "855547f13b6f4d901608c0c65e47d208393c8ced5baa0a99b1e625c939f49a41", + "b95c682aa102e852127cdd43f58ce4e381cf41e226b2698890eb279435f144ce", + "2257d311e774b07c5d50fa28a17faa8f50ab37fc493c54356fa3c5bf6a4245c4", + "5d4096a0ea1354d99d7847d6842791aca7fa9094a310e9b669a22c3765792c79", + "8a18d79b56ee00847942ff579ad9ce4a2e2220f33ebcbdd7222f0bde30c602f1", + "d34aa65128bff2fa149166b44deadf0175287558d33244548469f24c285cacb0" + ], + "nfldice": [ + "3a5bc1706db3ba2896eeb43c1e1d009b22000c05a4c9600c94a1820af264c213", + "2f5000122f933d2f259e3ea6f18da78939fe8025bec717d16b72edb709a7ccde", + "f9ac2fe4f0ca12ace763941fa743d6383a362e620a3d2d32ae8cef6871ca6d04", + "a2cd1f92e5e975a8ed18d1a4adb7c72435c0b4d957d753fc8f6b3ee0c6680dea", + "15f6c7e7775708ad7db5855d007f048e1ca1fb16e67107a6d87d31306d1d9846", + "d16d2caeb131cb8f8966bfcfc3ada9b6c6fa361a1f27182a0bfbe52e6395529f", + "841a3f2ed5102a2627ac9fee4c3036b702fe0d539db0290600032a7d8e57ea3d", + "e19316ed73de89865d8d476eb293798612be4c4c4aef2f9a1787454cd674a404", + "4a9ab14c74c76a606286398d549478c05f6c049a384e00b37dff908b0bb13e08", + "d8c325e8bb799f21cedcc03d6b7b464345741b46c30971fe7937bf69ecda4ed1" + ], + "betce": [ + "5b3c84885cb0a16834e3738c087feb80514fb27e8734dc473d8f572bed1a0aa3", + "b3520cac407d92db0cd44dea82398dcd3b3f0325f366f0893ad6478ee26a868a", + "2a7f5d74583aa2b9207bb435b414f70349d4a8769088f45cbd1a78df0ec22336", + "dfa7afe3aff927511e65d70007fe9af56be38b121d66982d9c7373c5cfcc5c19", + "4307374e19246ac522a6dfa10776cf78bbd9f9398bbdf255160da0299b9cf755", + "00089c17a13fc80bb98f3b4e75869330c9f38e9535493fa2ebe1be08b9c0eaa4", + "0733fa4d256e9edc6f3b62f914effba611034788bdb54f3af9d49da066906242", + "7ff46384c45a8569e36d4bf237effe9944b43c250dc3443d7ebe1a72010c07a1", + "de3db46c493fabdee83d4f079ede2169a0a900a13dce215c1923cfc5f1eb80c8", + "0fc3d29d5dca7d6226e6e7d6937d29ca61dd1161514d5792a5d2f3387ece82e0" + ], + "sfcef": [ + "68abcd48cfa8188ddf1d64e3afbcc6176bac8b5d0557e2512f345879fab208a7", + "613c9b53143362c7e81c09f52fa261ca0dd99e431f66be5960880d9f903ef7cf", + "62ff990b48f7249719d0cbc52fba5e2a28524889e152898f30c033793d70cc0c", + "2a7e2eb64e77e790bd14b5ea0225428780ffa40fecb4a8d5855f709dd8f07e9c", + "7ae38f0be0d2c1ba8d1ca617b931c0b1d324a0d8be33b6eff9fd3066235f4a5d", + "a9d7545c13d0a6f01f2c1bc665f7fbbfcfbc01a2bc01feb2078ccdc8861947cb", + "cb96df3b1a598df8f4a65f6f64db867eb3bd351840bafe25267d3ebac8cabcda", + "c82bfdea6c195a4a3f2fda07a50328feeb5d738378bfd12e353c0a6ae7f68868", + "3732962a2a54c7df89955780eb826e2deef4c998ef54cb1e939ce13bf65165cb", + "275de125f35a4aa1c90730afa74552ff3db4a39536b83651f7d657cc430004b8" + ] + }, + "sfcef_synthetic_n_diffs": { + "sfcef_synth_diag": 59, + "sfcef_synth_horiz": 58, + "sfcef_synth_vert": 29 + } +} diff --git a/tests/assets/sfcef_matlab_reference/sfcef_synth_diag_matlab.png b/tests/assets/sfcef_matlab_reference/sfcef_synth_diag_matlab.png new file mode 100644 index 0000000..1502bba Binary files /dev/null and b/tests/assets/sfcef_matlab_reference/sfcef_synth_diag_matlab.png differ diff --git a/tests/assets/sfcef_matlab_reference/sfcef_synth_horiz_matlab.png b/tests/assets/sfcef_matlab_reference/sfcef_synth_horiz_matlab.png new file mode 100644 index 0000000..c133b1b Binary files /dev/null and b/tests/assets/sfcef_matlab_reference/sfcef_synth_horiz_matlab.png differ diff --git a/tests/assets/sfcef_matlab_reference/sfcef_synth_vert_matlab.png b/tests/assets/sfcef_matlab_reference/sfcef_synth_vert_matlab.png new file mode 100644 index 0000000..c13707a Binary files /dev/null and b/tests/assets/sfcef_matlab_reference/sfcef_synth_vert_matlab.png differ diff --git a/tests/tools/utils.py b/tests/tools/utils.py index ad3aa55..6c0a661 100644 --- a/tests/tools/utils.py +++ b/tests/tools/utils.py @@ -9,7 +9,6 @@ import sqlite3 from typing import Any, Dict, Iterable, List, Tuple, Union, Literal, Optional, get_args import pickle -import copy import numpy as np from PIL import Image from shinier import utils @@ -20,7 +19,7 @@ # Constants & helpers # --------------------------------------------------------------------------- -# Mapping for utils.rgb2gray's conversion_type (when precomputing targets) +# Mapping for utils.rgb2gray's weighting_standard (when precomputing targets) AS_GRAY_NAME = {0: None, 1: "equal", 2: "rec601", 3: "rec709", 4: "rec2020"} # Get Image path @@ -43,6 +42,10 @@ ANSI_ESCAPE = re.compile(r'\x1B(?:[@-Z\\-_]|\[[0-?]*[ -/]*[@-~])') DB_PATH = Path(__file__).resolve().parent.parent / "hash_registry.db" +# Only terminal-success statuses are skipped on a re-run. 'failed'/'error' combos +# are intentionally NOT skipped so a real regression stays red until it is fixed. +TERMINAL_SKIP_STATUSES = ("done", "invalid", "invalid_option_combination") + # --------------------------------------------------------------------- # SQLite-backed HASH registry @@ -70,36 +73,50 @@ def _ensure_db() -> None: conn.commit() -def mark_hash_range_done(start: int, end: int) -> None: - """Mark all hashes whose integer value is in [start, end] as 'done'. +def mark_hash_range_done(start: int, end: int, namespace: str = "", chunk_size: int = 100_000) -> None: + """Upsert all integer hashes in [start, end] as 'done' (for START_AT skipping). Args: start: Lowest integer hash value (inclusive). end: Highest integer hash value (inclusive). + namespace: Prefix used when building hash keys (e.g. COVERAGE_MODE). + chunk_size: Rows per executemany batch. """ _ensure_db() if end < start: raise ValueError(f"end ({end}) must be >= start ({start})") + prefix = f"{namespace}:" if namespace else "" with sqlite3.connect(DB_PATH, timeout=30) as conn: - conn.execute( - """ - UPDATE hashes - SET status = 'done', - error = NULL, - timestamp = CURRENT_TIMESTAMP - WHERE CAST(hash AS INTEGER) BETWEEN ? AND ? - """, - (start, end), - ) + for batch_start in range(start, end + 1, chunk_size): + batch_end = min(end + 1, batch_start + chunk_size) + batch = [ + (f"{prefix}{i}", 'done') + for i in range(batch_start, batch_end) + ] + conn.executemany( + """ + INSERT INTO hashes (hash, status) VALUES (?, ?) + ON CONFLICT(hash) DO UPDATE SET + status = excluded.status, + timestamp = CURRENT_TIMESTAMP + """, + batch, + ) conn.commit() def is_already_done(hash_str: str) -> bool: - """Return True if the given combo hash is already processed (status is not 'pending').""" + """Return True only if the combo reached a terminal-success status. + + Combos marked 'failed' or 'error' are deliberately reported as *not* done so + that re-running the suite retries them; otherwise a single recorded failure + would be skipped forever and the test would go green without a fix. + """ _ensure_db() + placeholders = ",".join("?" * len(TERMINAL_SKIP_STATUSES)) with sqlite3.connect(DB_PATH) as conn: cur = conn.execute( - "SELECT 1 FROM hashes WHERE hash = ? AND status != 'pending' LIMIT 1", - (hash_str,), + f"SELECT 1 FROM hashes WHERE hash = ? AND status IN ({placeholders}) LIMIT 1", + (hash_str, *TERMINAL_SKIP_STATUSES), ) return cur.fetchone() is not None @@ -169,48 +186,10 @@ def get_all_hashes(status: Optional[str] = None) -> List[Tuple[str, str, str, Op return cur.fetchall() -def ensure_sequential_hashes( - max_hash: int, - *, - chunk_size: int = 100_000, -) -> None: - """Ensure the hash registry contains rows for '0'..str(max_hash). - - This function is idempotent: it uses INSERT OR IGNORE so re-running it - will only add missing rows and leave existing ones untouched. - - Args: - max_hash: Highest integer hash value to ensure, inclusive. The function - will ensure that every hash from "0" to str(max_hash) exists. - chunk_size: Number of rows to batch per executemany() call to avoid - building a gigantic parameter list in memory. - """ - _ensure_db() - with sqlite3.connect(DB_PATH, timeout=30) as conn: - cursor = conn.cursor() - for start in range(0, max_hash + 1, chunk_size): - end = min(max_hash + 1, start + chunk_size) - batch = [(str(i),) for i in range(start, end)] - cursor.executemany( - "INSERT OR IGNORE INTO hashes (hash) VALUES (?)", - batch, - ) - conn.commit() - -def initialize_db( - total_number_of_tests: int, - *, - chunk_size: int = 100_000, -) -> None: - """Drop and rebuild the hash registry with hashes 0..max_hash. - - Args: - total_number_of_tests: Highest integer hash value to create, inclusive. - chunk_size: Batch size for inserts. - """ +def initialize_db() -> None: + """Drop and recreate the hash registry (empty — rows are upserted on demand).""" reset_hash_registry(confirm=True) - ensure_sequential_hashes(max_hash=total_number_of_tests, chunk_size=chunk_size) def mark_hash_status( @@ -218,27 +197,24 @@ def mark_hash_status( status: str = "done", error: Optional[str] = None, ) -> None: - """Mark an existing hash row as 'done'. - - This is a lighter-weight helper than `register_hash` for the case where - the row is guaranteed to exist (e.g., after `ensure_sequential_hashes`). + """Upsert a hash row with the given status. Args: - hash_str: Hash identifier (e.g., "0", "1", "1180676"). - status: Status of the hash row to mark as 'done', 'invalide', etc. + hash_str: Hash identifier (e.g., "pruned:1180676"). + status: Status to record ('done', 'failed', 'invalid', etc.). error: Optional error message to record; defaults to NULL. """ _ensure_db() with sqlite3.connect(DB_PATH, timeout=30) as conn: conn.execute( """ - UPDATE hashes - SET status = ?, - error = ?, + INSERT INTO hashes (hash, status, error) VALUES (?, ?, ?) + ON CONFLICT(hash) DO UPDATE SET + status = excluded.status, + error = excluded.error, timestamp = CURRENT_TIMESTAMP - WHERE hash = ? """, - (status, error, hash_str), + (hash_str, status, error), ) conn.commit() @@ -500,34 +476,32 @@ def prepare_images(path_img: Path) -> Dict[str, Any]: Dict[str, Any]: Dictionary of image list (ImageListIO). """ out = {"buffers": {0: {0: {}, 1: {0: {}, 1: {}}}, 1: {0: {}, 1: {}}}} - out['buffers_other'] = copy.deepcopy(out['buffers']) out['images'] = None rec_standards = [r for r in get_args(REC_STANDARD)] images = ImageListIO(input_data=path_img) out['images'] = images - list_buffer = [np.zeros(im.shape, dtype=bool) for im in images] + # Memory: only the first image's color-treated buffer is ever consumed downstream + # (precompute_targets reads [...][0]). Building the full color-treated set for all 12 + # (as_gray, linear_luminance, rec) cells — plus the never-read buffers_other — would + # hold ~24 full-dataset copies (O(n_images * H * W), prohibitive at large resolutions). + # Treat a single-image collection per cell instead, so the footprint stays ~O(1 image). + first_mask = np.zeros(images[0].shape, dtype=bool) for ag in (0, 1): for ct in (0, 1): for rs in (1, 2, 3): - # ag, ct, rs = 0, 1, 2 - - # Prepare the images - buffers = ImageListIO(input_data=copy.deepcopy(list_buffer), conserve_memory=False) - buffers_other = ImageListIO(input_data=copy.deepcopy(list_buffer), conserve_memory=False) - for idx, image in enumerate(images): - buffers[idx] = image.astype(float) + buffers = ImageListIO(input_data=[first_mask.copy()], conserve_memory=False) + buffers[0] = images[0].astype(float) buffers.drange = (0, 255) - - # Convert them into xyY color space - rec_stardard = rec_standards[rs - 1] - output = ColorTreatment.forward_color_treatment(rec_standard=rec_stardard, input_images=buffers, - output_images=buffers, linear_luminance=bool(ct), - as_gray=ag, output_other=buffers_other) - _buffers, _buffers_other = output if isinstance(output, tuple) else (output, None) + # forward_color_treatment requires a sink for chroma channels on the + # color-preserving paths; we provide a throwaway and discard it. + other = ImageListIO(input_data=[first_mask.copy()], conserve_memory=False) + output = ColorTreatment.forward_color_treatment( + rec_standard=rec_standards[rs - 1], input_images=buffers, output_images=buffers, + linear_luminance=bool(ct), as_gray=ag, output_other=other) + _buffers = output[0] if isinstance(output, tuple) else output out["buffers"][ag][ct][rs] = _buffers - out["buffers_other"][ag][ct][rs] = _buffers_other return out diff --git a/tests/unit_tests/Converter_test.py b/tests/unit_tests/Converter_test.py index 5d0bcf7..e7f6e39 100644 --- a/tests/unit_tests/Converter_test.py +++ b/tests/unit_tests/Converter_test.py @@ -13,7 +13,8 @@ import numpy as np import pytest -from shinier.color import ColorConverter, WHITE_D65, COLOR_STANDARDS +from shinier.ImageListIO import ImageListIO +from shinier.color import ColorConverter, ColorTreatment, WHITE_D65, COLOR_STANDARDS, rgb2gray, rgb2ntsc_intensity pytestmark = pytest.mark.unit_tests @@ -140,4 +141,87 @@ def test_repr_and_assignment_behavior(): c = ColorConverter(rec_standard="rec709") c.rec_standard = "rec601" assert np.isclose(c.gamma, COLOR_STANDARDS["rec601"]["gamma"]) - assert np.allclose(c.M_RGB2XYZ, COLOR_STANDARDS["rec601"]["M_RGB2XYZ"]) \ No newline at end of file + assert np.allclose(c.M_RGB2XYZ, COLOR_STANDARDS["rec601"]["M_RGB2XYZ"]) + + +def test_rgb2ntsc_intensity_uses_matlab_yiq_weights(rgb_sample: np.ndarray): + """NTSC intensity is the Y channel from MATLAB's rgb2ntsc/YIQ transform.""" + weights = np.array([0.298936021293775, 0.587043074451121, 0.114020904255103]) + direct = rgb2gray(rgb_sample, weighting_standard="rec601", matlab_601=True) + + np.testing.assert_allclose(direct, np.tensordot(rgb_sample, weights, axes=([-1], [0]))) + np.testing.assert_allclose(rgb2ntsc_intensity(rgb_sample), direct) + + +def test_legacy_grayscale_forward_uses_matlab_ntsc_intensity(rgb_sample: np.ndarray): + """Legacy grayscale color treatment uses MATLAB-compatible NTSC/YIQ intensity.""" + image = (rgb_sample * 255).astype(np.float64) + images = ImageListIO(input_data=[image], conserve_memory=False) + + result, other = ColorTreatment.forward_color_treatment( + rec_standard="rec601", + input_images=images, + output_images=images, + linear_luminance=False, + as_gray=True, + legacy_mode=True, + ) + + np.testing.assert_allclose(result[0], rgb2ntsc_intensity(image)) + assert other is None + + +def test_legacy_grayscale_backward_does_not_gamma_encode() -> None: + """Legacy grayscale output keeps MATLAB-compatible intensities unchanged.""" + gray = np.array([[0.0, 64.0], [128.0, 255.0]], dtype=np.float64) + images = ImageListIO(input_data=[gray], conserve_memory=False) + + result = ColorTreatment.backward_color_treatment( + rec_standard="rec601", + input_images=images, + output_images=images, + linear_luminance=False, + as_gray=True, + legacy_mode=True, + ) + + np.testing.assert_allclose(result[0], np.dstack([gray, gray, gray])) + + +def test_legacy_grayscale_backward_skips_gamma_non_legacy_applies_it() -> None: + """Explicitly verify that legacy_mode controls gamma encoding in the grayscale backward pass. + + legacy_mode=True must return Y values unchanged (no sRGB transfer function). + legacy_mode=False must apply linRGB_to_sRGB, producing different values for mid-range inputs. + If the legacy_mode branch is missing or broken both results would be identical. + """ + from shinier.color.Converter import ColorConverter + # Mid-range linear values where gamma makes a measurable difference (not 0 or 1) + gray = np.array([[50.0, 100.0], [150.0, 200.0]], dtype=np.float64) + + images_legacy = ImageListIO(input_data=[gray.copy()], conserve_memory=False) + images_non_legacy = ImageListIO(input_data=[gray.copy()], conserve_memory=False) + + result_legacy = ColorTreatment.backward_color_treatment( + rec_standard="rec601", input_images=images_legacy, output_images=images_legacy, + linear_luminance=False, as_gray=True, legacy_mode=True, + ) + result_non_legacy = ColorTreatment.backward_color_treatment( + rec_standard="rec601", input_images=images_non_legacy, output_images=images_non_legacy, + linear_luminance=False, as_gray=True, legacy_mode=False, + ) + + # legacy_mode=True: output must equal input replicated (no gamma applied) + np.testing.assert_allclose(result_legacy[0], np.dstack([gray, gray, gray]), atol=1e-6, + err_msg="legacy_mode=True must leave Y values unchanged") + + # legacy_mode=False: output must equal gamma-encoded Y replicated + converter = ColorConverter(rec_standard="rec601") + Y_linear = gray / 255.0 + Y_gamma = converter.linRGB_to_sRGB(Y_linear[..., np.newaxis])[..., 0] * 255 + np.testing.assert_allclose(result_non_legacy[0], np.dstack([Y_gamma, Y_gamma, Y_gamma]), atol=1e-6, + err_msg="legacy_mode=False must apply sRGB gamma encoding") + + # The two must differ — if the branch was missing this would fail + assert not np.allclose(result_legacy[0], result_non_legacy[0]), \ + "legacy_mode=True and legacy_mode=False must produce different outputs for mid-range values" diff --git a/tests/unit_tests/ImageListIO_test.py b/tests/unit_tests/ImageListIO_test.py index 360c61d..b94ae2d 100644 --- a/tests/unit_tests/ImageListIO_test.py +++ b/tests/unit_tests/ImageListIO_test.py @@ -66,11 +66,12 @@ def test_as_gray_flag_validation() -> None: ImageListIO(input_data=arrays, as_gray=9) # invalid literal -def test_as_gray_flag_accepted(test_tmpdir: Path) -> None: +@pytest.mark.parametrize("as_gray", [3, 4]) +def test_as_gray_flag_accepted(test_tmpdir: Path, as_gray: int) -> None: """as_gray flag should be correctly forwarded and retained.""" arrays = [_make_rgb(seed=s) for s in range(2)] - coll = ImageListIO(input_data=arrays, conserve_memory=True, as_gray=3) - assert coll.as_gray == 3 + coll = ImageListIO(input_data=arrays, conserve_memory=True, as_gray=as_gray) + assert coll.as_gray == as_gray coll.close() @@ -140,4 +141,4 @@ def test_readonly_copy_behaves_identically(test_tmpdir: Path) -> None: with pytest.raises(RuntimeError): clone[0] = _make_rgb() # should not allow modification coll.close() - clone.close() \ No newline at end of file + clone.close() diff --git a/tests/unit_tests/ImageProcessor_test.py b/tests/unit_tests/ImageProcessor_test.py index 9f703ed..b0981d0 100644 --- a/tests/unit_tests/ImageProcessor_test.py +++ b/tests/unit_tests/ImageProcessor_test.py @@ -105,6 +105,36 @@ def test_processor_init_from_folder(test_tmpdir: Path): assert len(proc.dataset.images) > 0 +def test_print_log_results_records_options(test_tmpdir: Path) -> None: + """Processing logs should include the selected options.""" + inp = _prepare_temp_images(test_tmpdir) + out = test_tmpdir / "OUTPUT" + out.mkdir() + opt = Options(input_folder=inp, output_folder=out, mode=9, standalone_op="ie_methods", ie_methods="tidhe") + proc = ImageProcessor(dataset=ImageDataset(options=opt), options=opt, verbose=0, from_unit_test=True) + proc.validation.append({ + "iter": 0, + "step": 0, + "processing_function": "hist_match", + "image": "test.png", + "channel": None, + "valid_result": True, + "log_result": "RMS error = 0", + }) + + proc.print_log_results() + + logs = sorted(out.glob("log_*.txt")) + assert len(logs) == 1 + text = logs[0].read_text() + assert "[Options]" in text + assert "mode: 9" in text + assert "standalone_op: 'ie_methods'" in text + assert "ie_methods: 'tidhe'" in text + assert f"output_folder: {out}" in text + assert "channel=all; processing function=hist_match" in text + + def test_processor_log_and_results(test_tmpdir: Path) -> None: """Processor should produce log entries and valid result arrays.""" arrays = [_make_rgb(seed=s) for s in range(2)] @@ -365,3 +395,155 @@ def test_lum_match_constant_image_zero_std_is_stable(test_tmpdir: Path) -> None: mean_after, std_after = _stats(proc.dataset.buffer[0]) assert mean_after == pytest.approx(140.0, abs=1e-6) assert std_after == pytest.approx(0.0, abs=1e-9) + + +def _build_hist_match_processor( + test_tmpdir: Path, + images: list[np.ndarray], + target_hist: str, + as_gray: bool = True, +) -> ImageProcessor: + """Build a processor wired for direct hist_match unit tests.""" + inp, out = _make_input_output_dirs(test_tmpdir) + opt = Options( + input_folder=inp, + output_folder=out, + mode=2, + target_hist=target_hist, + as_gray=as_gray, + hist_optim=False, + verbose=-1, + ) + io = ImageListIO(input_data=images, conserve_memory=True) + ds = ImageDataset(images=io, options=opt) + ds.buffer = ImageListIO( + input_data=[im.astype(np.float64) for im in images], conserve_memory=True + ) + proc = ImageProcessor(dataset=ds, options=opt, verbose=-1, from_unit_test=True) + proc.bool_masks = [np.ones(im.shape[:2], dtype=bool) for im in images] + proc._sum_bool_masks = [ + np.ones((*im.shape[:2], 1), dtype=bool) for im in images + ] + proc._is_first_operation = True + proc._is_last_operation = True + proc._initial_buffer = ds.buffer + proc._compute_initial_target_histogram() + return proc + + +def _build_ie_processor( + test_tmpdir: Path, + images: list[np.ndarray], + ie_algorithm: str = "classic_he", +) -> ImageProcessor: + """Build a processor wired for direct ie_methods unit tests (mode=9).""" + inp, out = _make_input_output_dirs(test_tmpdir) + opt = Options( + input_folder=inp, + output_folder=out, + mode=9, + standalone_op="ie_methods", + ie_methods=ie_algorithm, + as_gray=True, + verbose=-1, + ) + io = ImageListIO(input_data=images, conserve_memory=True) + ds = ImageDataset(images=io, options=opt) + ds.buffer = ImageListIO( + input_data=[im.astype(np.float64) for im in images], conserve_memory=True + ) + proc = ImageProcessor(dataset=ds, options=opt, verbose=-1, from_unit_test=True) + proc.bool_masks = [np.ones(im.shape[:2], dtype=bool) for im in images] + proc._sum_bool_masks = [ + np.ones((*im.shape[:2], 1), dtype=bool) for im in images + ] + proc._is_first_operation = True + proc._is_last_operation = True + proc._initial_buffer = ds.buffer + return proc + + +def _make_gray(h: int = 64, w: int = 64, seed: int = 0) -> np.ndarray: + """Create a random grayscale image with float64 values in [0, 255].""" + rng = np.random.default_rng(seed) + return rng.integers(0, 256, (h, w), dtype=np.uint8).astype(np.float64) + + +@pytest.mark.parametrize("ie_algorithm", ["classic_he", "tidhe", "rdfhe", "nfldice", "betce", "sfcef"]) +def test_image_enhancement_output_shape_preserved(test_tmpdir: Path, ie_algorithm: str) -> None: + """ie_methods must preserve image shape for each image-enhancement algorithm.""" + img3d = _make_gray()[:, :, np.newaxis] + proc = _build_ie_processor(test_tmpdir, [img3d], ie_algorithm=ie_algorithm) + proc.ie_methods() + assert proc.dataset.buffer[0].shape == img3d.shape + + +@pytest.mark.parametrize("ie_algorithm", ["classic_he", "tidhe", "rdfhe", "nfldice", "betce", "sfcef"]) +def test_image_enhancement_output_range(test_tmpdir: Path, ie_algorithm: str) -> None: + """ie_methods output must stay within [0, 255].""" + img3d = _make_gray()[:, :, np.newaxis] + proc = _build_ie_processor(test_tmpdir, [img3d], ie_algorithm=ie_algorithm) + proc.ie_methods() + result = proc.dataset.buffer[0] + assert result.min() >= -1e-6 + assert result.max() <= 255 + 1e-6 + + +def test_image_enhancement_records_article_metrics(test_tmpdir: Path) -> None: + """ie_methods should log diagnostic metrics used to evaluate enhancement.""" + img3d = _make_gray()[:, :, np.newaxis] + proc = _build_ie_processor(test_tmpdir, [img3d], ie_algorithm="tidhe") + + proc.ie_methods() + + log_result = proc.validation[-1]["log_result"] + assert proc.validation[-1]["processing_function"] == "ie_methods" + assert proc.validation[-1]["channel"] == "all" + for metric in ("AMBE=", "MSSIM=", "PSNR=", "BP2BPSIM=", "CI=", "Entropy="): + assert metric in log_result + + +@pytest.mark.parametrize("ie_algorithm", ["classic_he", "tidhe", "rdfhe", "nfldice", "betce", "sfcef"]) +def test_image_enhancement_multiple_images_independent(test_tmpdir: Path, ie_algorithm: str) -> None: + """ie_methods must process each image independently.""" + img1 = np.full((64, 64, 1), 80.0, dtype=np.float64) + img1[10:30, 10:30] = 200.0 + img2 = _make_gray(seed=7)[:, :, np.newaxis] + proc = _build_ie_processor(test_tmpdir, [img1, img2], ie_algorithm=ie_algorithm) + proc.ie_methods() + for r in (proc.dataset.buffer[0], proc.dataset.buffer[1]): + assert r.min() >= -1e-6 + assert r.max() <= 255 + 1e-6 + + +@pytest.mark.parametrize("ie_algorithm", ["classic_he", "tidhe", "rdfhe", "nfldice", "betce", "sfcef"]) +def test_image_enhancement_runs_with_mask_present(test_tmpdir: Path, ie_algorithm: str) -> None: + """ie_methods must still run when masks are attached to the processor.""" + img = _make_gray()[:, :, np.newaxis] + proc = _build_ie_processor(test_tmpdir, [img], ie_algorithm=ie_algorithm) + mask = np.zeros((64, 64), dtype=bool) + mask[:32, :] = True + proc.bool_masks = [mask] + proc.ie_methods() + result = proc.dataset.buffer[0] + assert result.shape == img.shape + assert result.min() >= -1e-6 + assert result.max() <= 255 + 1e-6 + + +def test_equal_target_hist_accepted_and_runs(test_tmpdir: Path) -> None: + """target_hist='equal' must be accepted by Options and produce a valid result.""" + img = _make_gray()[:, :, np.newaxis] + proc = _build_hist_match_processor(test_tmpdir, [img], "equal") + proc.hist_match() + result = proc.dataset.buffer[0] + assert result.shape == img.shape + assert result.min() >= -1e-6 + assert result.max() <= 255 + 1e-6 + + +def test_target_hist_unrecognised_string_rejected_by_options(test_tmpdir: Path) -> None: + """Options must reject any unrecognised target_hist string.""" + inp, out = _make_input_output_dirs(test_tmpdir) + with pytest.raises((ValueError, ValidationError)): + Options(input_folder=inp, output_folder=out, target_hist="not_a_valid_value") diff --git a/tests/unit_tests/Options_test.py b/tests/unit_tests/Options_test.py index 5d09b66..2866542 100644 --- a/tests/unit_tests/Options_test.py +++ b/tests/unit_tests/Options_test.py @@ -220,6 +220,53 @@ def test_target_spectrum_validation(tmp_dirs): Options(input_folder=in_dir, output_folder=out_dir, target_spectrum=bad_path_type) +@pytest.mark.unit_tests +def test_target_hist_equal_accepted(tmp_dirs): + """target_hist='equal' must be accepted and stored as-is.""" + in_dir, out_dir, _ = tmp_dirs + opt = Options(input_folder=in_dir, output_folder=out_dir, target_hist="equal") + assert opt.target_hist == "equal" + + +@pytest.mark.unit_tests +@pytest.mark.parametrize("ie_methods", ["classic_he", "tidhe", "rdfhe", "nfldice", "betce", "sfcef"]) +def test_ie_methods_accepted(tmp_dirs, ie_methods): + """All image-enhancement methods must be accepted in mode 9.""" + in_dir, out_dir, _ = tmp_dirs + opt = Options(input_folder=in_dir, output_folder=out_dir, + mode=9, standalone_op="ie_methods", ie_methods=ie_methods) + assert opt.standalone_op == "ie_methods" + assert opt.ie_methods == ie_methods + + +@pytest.mark.unit_tests +def test_operation_dithering_requires_dithering_nonzero(tmp_dirs): + """mode=9 with standalone_op='dithering' and dithering=0 must raise.""" + in_dir, out_dir, _ = tmp_dirs + with pytest.raises((ValueError, ValidationError)): + Options(input_folder=in_dir, output_folder=out_dir, + mode=9, standalone_op="dithering", dithering=0) + + +@pytest.mark.unit_tests +def test_mode9_ie_methods_rejects_dithering(tmp_dirs): + """mode=9 image-enhancement methods must not accept dithering.""" + in_dir, out_dir, _ = tmp_dirs + with pytest.raises((ValueError, ValidationError)): + Options(input_folder=in_dir, output_folder=out_dir, + mode=9, standalone_op="ie_methods", dithering=2) + + +@pytest.mark.unit_tests +def test_target_hist_unrecognised_string_rejected(tmp_dirs): + """Arbitrary strings that are not valid sentinel values must be rejected.""" + in_dir, out_dir, _ = tmp_dirs + with pytest.raises((ValueError, ValidationError)): + Options(input_folder=in_dir, output_folder=out_dir, target_hist="flat") + with pytest.raises((ValueError, ValidationError)): + Options(input_folder=in_dir, output_folder=out_dir, target_hist="equalization") + + @pytest.mark.unit_tests def test_hist_optim_overwrites_hist_spec(tmp_dirs): """hist_optim=True should nullify hist_specification.""" @@ -239,10 +286,11 @@ def test_rescaling_forbidden_modes_overwrite(tmp_dirs): @pytest.mark.unit_tests def test_mode9_dithering_zero_raises(tmp_dirs): - """Mode 9 cannot have dithering=0.""" + """Mode 9 with standalone_op='dithering' and dithering=0 must raise.""" in_dir, out_dir, _ = tmp_dirs - with pytest.raises(ValueError): - Options(input_folder=in_dir, output_folder=out_dir, mode=9, dithering=0) + with pytest.raises((ValueError, ValidationError)): + Options(input_folder=in_dir, output_folder=out_dir, + mode=9, standalone_op="dithering", dithering=0) @pytest.mark.unit_tests @@ -271,6 +319,8 @@ def test_legacy_mode_overrides(tmp_dirs): opt = Options(input_folder=in_dir, output_folder=out_dir, legacy_mode=True) assert not opt.conserve_memory assert opt.as_gray + assert not opt.linear_luminance + assert opt.rec_standard == 1 assert opt.dithering == 0 assert opt.hist_specification == 1 assert not opt.safe_lum_match @@ -340,8 +390,11 @@ def test_all_combo(tmp_dirs): continue params = dict(zip(keys, combo)) - if params['mode'] == 9 and params['dithering'] == 0: - params['dithering'] = 1 + if params['mode'] == 9: + if params['standalone_op'] == 'dithering' and params['dithering'] == 0: + params['dithering'] = 1 + elif params['standalone_op'] == 'ie_methods' and params['dithering'] != 0: + params['dithering'] = 0 Options(**params) if pbar is not None: pbar.update(1) diff --git a/tests/unit_tests/PrunedCoverage_test.py b/tests/unit_tests/PrunedCoverage_test.py new file mode 100644 index 0000000..ebac1e7 --- /dev/null +++ b/tests/unit_tests/PrunedCoverage_test.py @@ -0,0 +1,109 @@ +"""Unit coverage for parameters collapsed in the pruned validation mode. + +Each test exercises one value that the pruned mode skips, verifying that it runs +without error and produces finite output. See tests/README.md for rationale. +""" +import numpy as np +import pytest +from pathlib import Path +from shinier import ImageDataset, Options, ImageProcessor +from tests import utils as utils_test + +pytestmark = pytest.mark.unit_tests + +_BASE_OPTS = dict( + mode=2, + as_gray=False, + linear_luminance=False, + rec_standard=2, + legacy_mode=False, + iterations=1, + hist_optim=False, + hist_specification=1, + rescaling=2, + dithering=0, + conserve_memory=False, + fft_padding_mode=0, + standalone_op="ie_methods", + ie_methods="tidhe", + seed=None, + verbose=-1, +) + + +def _run(tmp_path: Path, **overrides) -> ImageProcessor: + tmp_path.mkdir(parents=True, exist_ok=True) + buffers = utils_test.prepare_images(utils_test.IMAGE_PATH) + images = utils_test.select_n_imgs(buffers["images"], n=2, seed=0) + opts = Options(**{**_BASE_OPTS, "output_folder": tmp_path, **overrides}) + ds = ImageDataset.model_construct(images=images, options=opts) + return ImageProcessor.model_construct(dataset=ds, options=opts, verbose=-1, from_validation_test=True) + + +# ----- ie_methods (pruned keeps only "tidhe") ----- + +@pytest.mark.parametrize("method", ["classic_he", "rdfhe", "nfldice", "betce", "sfcef"]) +def test_ie_method_runs(method, tmp_path): + proc = _run(tmp_path, mode=9, standalone_op="ie_methods", ie_methods=method) + for img in proc._final_buffer: + assert np.all(np.isfinite(img)), f"{method}: NaN/Inf in output" + + +# ----- fft_padding_mode (pruned keeps 0 and 1) ----- + +@pytest.mark.parametrize("padding_mode", [2, 3]) +def test_fft_padding_mode_runs(padding_mode, tmp_path): + proc = _run(tmp_path, mode=3, fft_padding_mode=padding_mode) + for img in proc._final_buffer: + assert np.all(np.isfinite(img)), f"fft_padding_mode={padding_mode}: NaN/Inf in output" + + +# ----- rec_standard (pruned keeps only 1 = rec601) ----- + +@pytest.mark.parametrize("rs", [2, 3]) +def test_rec_standard_runs(rs, tmp_path): + proc = _run(tmp_path, rec_standard=rs) + assert proc._final_buffer is not None + + +# ----- dithering (pruned keeps 0 and 1) ----- + +def test_dithering_2_floyd_steinberg_runs(tmp_path): + proc = _run(tmp_path, dithering=2, as_gray=True) + for img in proc._final_buffer: + assert np.all(np.isfinite(img)) + + +# ----- hist_specification (pruned keeps 1 and 3) ----- + +@pytest.mark.parametrize("hs", [2, 4]) +def test_hist_specification_runs(hs, tmp_path): + proc = _run(tmp_path, mode=2, hist_specification=hs) + assert proc._final_buffer is not None + + +# ----- rescaling (pruned keeps 0 and 2) ----- + +@pytest.mark.parametrize("rescaling", [1, 3]) +def test_rescaling_runs(rescaling, tmp_path): + proc = _run(tmp_path, mode=3, rescaling=rescaling) + for img in proc._final_buffer: + assert np.all(np.isfinite(img)), f"rescaling={rescaling}: NaN/Inf in output" + + +# ----- conserve_memory (pruned keeps False) ----- + +def test_conserve_memory_true_matches_false(tmp_path): + proc_off = _run(tmp_path / "off", conserve_memory=False, seed=12345) + proc_on = _run(tmp_path / "on", conserve_memory=True, seed=12345) + for a, b in zip(proc_off._final_buffer, proc_on._final_buffer): + np.testing.assert_allclose(a, b, atol=1e-6) + + +# ----- seed (pruned keeps None) ----- + +def test_seeded_run_is_reproducible(tmp_path): + proc1 = _run(tmp_path / "r1", mode=2, hist_specification=2, seed=9999) + proc2 = _run(tmp_path / "r2", mode=2, hist_specification=2, seed=9999) + for a, b in zip(proc1._final_buffer, proc2._final_buffer): + np.testing.assert_array_equal(a, b) diff --git a/tests/unit_tests/SHINIER_CLI_test.py b/tests/unit_tests/SHINIER_CLI_test.py index d067783..8f81593 100644 --- a/tests/unit_tests/SHINIER_CLI_test.py +++ b/tests/unit_tests/SHINIER_CLI_test.py @@ -7,55 +7,143 @@ pytestmark = pytest.mark.unit_tests -def test_prompt_bool_retries_after_invalid_input(monkeypatch, capsys, recwarn): - """Invalid bool input should repeat the same prompt and preserve the default.""" - answers = iter(["¥", "n"]) +def _mock_input(monkeypatch, answers): + answers = iter(answers) monkeypatch.setattr(builtins, "input", lambda: next(answers)) - value = prompt("Safe luminance matching?", default="n", kind="bool") +def _assert_invalid_retries(monkeypatch, capsys, label, answers, expected, **kwargs): + _mock_input(monkeypatch, answers) + value = prompt(label, **kwargs) captured = capsys.readouterr() - assert value is False + + assert value == expected assert "Invalid input." in captured.out - assert captured.out.count("Safe luminance matching?") == 2 + assert captured.out.count(label) == 2 + + +@pytest.mark.parametrize( + "default, expected", + [("n", False), ("y", True)], +) +def test_prompt_bool_default_returns_boolean(monkeypatch, default, expected): + """Pressing Enter for a bool prompt should return True/False.""" + _mock_input(monkeypatch, [""]) + + assert prompt("Safe luminance matching?", default=default, kind="bool") is expected + + +def test_prompt_bool_retries_after_invalid_input(monkeypatch, capsys, recwarn): + """Invalid bool input should re-prompt without warning.""" + _assert_invalid_retries( + monkeypatch, capsys, "Safe luminance matching?", ["¥", "n"], + False, default="n", kind="bool" + ) assert len(recwarn) == 0 -def test_prompt_bool_default_returns_boolean(monkeypatch): - """Pressing Enter for a bool prompt should return True/False, not a choice index.""" - monkeypatch.setattr(builtins, "input", lambda: "") +def test_prompt_choice_default_returns_index(monkeypatch): + """Choice defaults should keep their existing index-based return value.""" + _mock_input(monkeypatch, [""]) - value = prompt("Safe luminance matching?", default="n", kind="bool") + assert prompt("Processing mode", default=2, kind="choice", choices=["A", "B"]) == 2 - assert value is False +@pytest.mark.parametrize( + "answers, expected", + [(["bad", "2"], 2), (["99", "1"], 1)], +) +def test_prompt_choice_invalid_input_retries(monkeypatch, capsys, answers, expected): + """Invalid choice input should re-prompt, not quit.""" + _assert_invalid_retries( + monkeypatch, capsys, "Processing mode", answers, expected, + default=1, kind="choice", choices=["A", "B"] + ) -def test_prompt_bool_yes_default_returns_boolean(monkeypatch): - """Bool defaults should preserve boolean semantics for yes defaults too.""" - monkeypatch.setattr(builtins, "input", lambda: "") - value = prompt("Histogram optimization?", default="y", kind="bool") +def test_prompt_int_enter_returns_default(monkeypatch): + """Enter on an int prompt should return the default integer.""" + _mock_input(monkeypatch, [""]) - assert value is True + value = prompt("Kernel size", default=3, kind="int") + assert value == 3 + assert isinstance(value, int) -def test_prompt_choice_retries_after_invalid_input(monkeypatch, capsys): - """Invalid choice input should repeat the prompt with the same default.""" - answers = iter(["bad", "2"]) - monkeypatch.setattr(builtins, "input", lambda: next(answers)) - value = prompt("Processing mode", default=1, kind="choice", choices=["A", "B"]) +@pytest.mark.parametrize( + "answers, expected, kwargs", + [ + (["abc", "7"], 7, {}), + (["200", "5"], 5, {"min_v": 1, "max_v": 10}), + ], +) +def test_prompt_int_invalid_input_retries(monkeypatch, capsys, answers, expected, kwargs): + """Invalid int input should re-prompt, not quit.""" + _assert_invalid_retries( + monkeypatch, capsys, "Kernel size", answers, expected, + default=3, kind="int", **kwargs + ) - captured = capsys.readouterr() - assert value == 2 - assert "Invalid input." in captured.out - assert captured.out.count("Processing mode") == 2 +def test_prompt_float_enter_returns_default(monkeypatch): + """Enter on a float prompt should return the default float.""" + _mock_input(monkeypatch, [""]) -def test_prompt_choice_default_returns_index(monkeypatch): - """Choice defaults should keep their existing index-based return value.""" - monkeypatch.setattr(builtins, "input", lambda: "") + value = prompt("Learning rate", default=0.5, kind="float") - value = prompt("Processing mode", default=2, kind="choice", choices=["A", "B"]) + assert value == 0.5 + assert isinstance(value, float) + + +@pytest.mark.parametrize( + "label, answers, expected, kwargs", + [ + ("Learning rate", ["notanumber", "1.2"], 1.2, {}), + ("Clip threshold", ["-1.0", "0.8"], 0.8, {"min_v": 0.0, "max_v": 1.0}), + ], +) +def test_prompt_float_invalid_input_retries(monkeypatch, capsys, label, answers, expected, kwargs): + """Invalid float input should re-prompt, not quit.""" + _assert_invalid_retries( + monkeypatch, capsys, label, answers, pytest.approx(expected), + default=0.5, kind="float", **kwargs + ) + + +def test_prompt_str_enter_returns_default(monkeypatch): + """Enter on a str prompt should return the default string.""" + _mock_input(monkeypatch, [""]) + + assert prompt("Output name", default="result") == "result" + + +def test_prompt_tuple_enter_returns_default(monkeypatch): + """Enter on a tuple prompt should return the default tuple.""" + _mock_input(monkeypatch, [""]) + + assert prompt("Range", default=(0.0, 1.0), kind="tuple") == (0.0, 1.0) + + +def test_prompt_tuple_invalid_retries(monkeypatch, capsys): + """Malformed tuple input should re-prompt, not quit.""" + _assert_invalid_retries( + monkeypatch, capsys, "Range", ["notvalid,,", "0.0, 1.0"], + (0.0, 1.0), default=(0.0, 1.0), kind="tuple" + ) + + +def test_prompt_validator_failure_retries(monkeypatch, capsys): + """A failing validator should re-prompt instead of returning the bad value.""" + def validator(v): + if v == "bad_path": + return False, "Path does not exist" + return True, "" + + _mock_input(monkeypatch, ["bad_path", "good_path"]) + value = prompt("Input path", default="default_path", validator=validator) + captured = capsys.readouterr() - assert value == 2 + assert value == "good_path" + assert "Path does not exist" in captured.out + assert captured.out.count("Input path") == 2 diff --git a/tests/unit_tests/Utils_test.py b/tests/unit_tests/Utils_test.py new file mode 100644 index 0000000..23b5b9a --- /dev/null +++ b/tests/unit_tests/Utils_test.py @@ -0,0 +1,163 @@ +import numpy as np +import pytest + +from shinier.utils import ( + betce_gray, + classic_he_gray, + compute_ambe, + compute_bp2bpsim, + compute_contrast_improvement, + compute_image_entropy, + compute_mssim, + compute_psnr, + nfldice_gray, + rdfhe_gray, + sfcef_gray, + tidhe_gray, +) + +pytestmark = pytest.mark.unit_tests + +IE_METHODS = [ + classic_he_gray, + tidhe_gray, + rdfhe_gray, + nfldice_gray, + betce_gray, + sfcef_gray, +] +LUT_METHODS = [classic_he_gray, tidhe_gray, rdfhe_gray, nfldice_gray, betce_gray] + + +@pytest.mark.parametrize("method", IE_METHODS) +def test_image_enhancement_rejects_non_grayscale(method) -> None: + """Image-enhancement methods are intentionally limited to 2D grayscale images.""" + with pytest.raises(ValueError, match="2D grayscale"): + method(np.zeros((4, 4, 1), dtype=np.uint8)) + + +@pytest.mark.parametrize("method", IE_METHODS) +def test_image_enhancement_rejects_non_uint8(method) -> None: + """Image-enhancement methods use the MATLAB uint8 reference domain.""" + with pytest.raises(ValueError, match="dtype uint8"): + method(np.zeros((4, 4), dtype=np.float64)) + + +@pytest.mark.parametrize("method", IE_METHODS) +def test_image_enhancement_output_shape_dtype_and_range(method) -> None: + rng = np.random.default_rng(42) + image = rng.integers(0, 256, (64, 64), dtype=np.uint8) + result = method(image) + assert result.shape == image.shape + assert result.dtype == np.uint8 + assert np.min(result) >= 0 + assert np.max(result) <= 255 + + +@pytest.mark.parametrize("method", IE_METHODS) +def test_image_enhancement_constant_image_stable(method) -> None: + result = method(np.full((64, 64), 128, dtype=np.uint8)) + assert result.shape == (64, 64) + assert result.dtype == np.uint8 + assert np.min(result) >= 0 + assert np.max(result) <= 255 + + +@pytest.mark.parametrize("method", LUT_METHODS) +def test_image_enhancement_mapping_is_monotonic(method) -> None: + image = np.arange(256, dtype=np.uint8).reshape(16, 16) + result = method(image) + assert np.all(np.diff(result.ravel().astype(np.int64)) >= 0) + + +def test_classic_he_gray_constant_image_unchanged() -> None: + image = np.full((16, 16), 128, dtype=np.uint8) + np.testing.assert_array_equal(classic_he_gray(image), image) + + +def test_classic_he_gray_two_levels_use_full_range() -> None: + image = np.zeros((16, 16), dtype=np.uint8) + image[:, 8:] = 128 + result = classic_he_gray(image) + assert np.min(result) == 0 + assert np.max(result) == 255 + + +def test_nfldice_gray_midpoint_maps_to_half_range() -> None: + """With p_l=127.5, level 128 maps just above 0.5*(L-1), matching MATLAB.""" + image = np.full((8, 8), 128, dtype=np.uint8) + assert int(nfldice_gray(image)[0, 0]) == 129 + + +def test_contrast_metric_ambe_known_values() -> None: + image = np.zeros((4, 4), dtype=np.uint8) + shifted = np.full((4, 4), 10, dtype=np.uint8) + assert compute_ambe(image, image) == pytest.approx(0.0) + assert compute_ambe(image, shifted) == pytest.approx(10.0) + + +def test_contrast_metric_ambe_3d_averages_channels() -> None: + reference = np.zeros((2, 2, 2), dtype=np.float64) + enhanced = np.zeros_like(reference) + enhanced[:, :, 0] = 10 + enhanced[:, :, 1] = -10 + assert compute_ambe(reference, enhanced) == pytest.approx(10.0) + + +def test_contrast_metric_ci_and_entropy_constant_and_varied() -> None: + constant = np.zeros((4, 4), dtype=np.uint8) + varied = np.array([[0, 1], [0, 1]], dtype=np.uint8) + + assert compute_contrast_improvement(constant, n_bins=2) == pytest.approx(0.0) + assert compute_image_entropy(constant, n_bins=2) == pytest.approx(0.0) + assert compute_contrast_improvement(varied, n_bins=2) > 0 + assert compute_image_entropy(varied, n_bins=2) == pytest.approx(1.0) + + +def test_contrast_metric_mssim_identical_is_one() -> None: + image = np.arange(64, dtype=np.float64).reshape(8, 8) + assert compute_mssim(image, image, data_range=63) == pytest.approx(1.0) + + +def test_contrast_metric_psnr_known_values() -> None: + reference = np.zeros((4, 4), dtype=np.float64) + enhanced = np.ones((4, 4), dtype=np.float64) + + assert compute_psnr(reference, reference, data_range=1) == np.inf + assert compute_psnr(reference, enhanced, data_range=1) == pytest.approx(0.0) + + +def test_contrast_metric_psnr_3d_averages_channels() -> None: + reference = np.zeros((2, 2, 2), dtype=np.float64) + enhanced = np.zeros_like(reference) + enhanced[:, :, 0] = 1 + enhanced[:, :, 1] = 2 + expected = (0.0 + (-10.0 * np.log10(4.0))) / 2.0 + assert compute_psnr(reference, enhanced, data_range=1) == pytest.approx(expected) + + +def test_contrast_metric_bp2bpsim_known_values() -> None: + reference = np.array([[0]], dtype=np.uint8) + enhanced = np.array([[1]], dtype=np.uint8) + + assert compute_bp2bpsim(reference, reference) == pytest.approx(1.0) + assert compute_bp2bpsim(reference, enhanced) == pytest.approx(7.0 / 8.0) + assert compute_bp2bpsim(reference, enhanced, n_bits=1) == pytest.approx(0.0) + + +def test_contrast_metric_bp2bpsim_3d_averages_channels() -> None: + reference = np.zeros((1, 1, 2), dtype=np.uint8) + enhanced = np.zeros_like(reference) + enhanced[:, :, 0] = 1 + assert compute_bp2bpsim(reference, enhanced, n_bits=1) == pytest.approx(0.5) + + +@pytest.mark.parametrize("metric", [compute_ambe, compute_mssim, compute_psnr, compute_bp2bpsim]) +def test_pairwise_contrast_metrics_reject_shape_mismatch(metric) -> None: + with pytest.raises(ValueError, match="same shape"): + metric(np.zeros((4, 4)), np.zeros((4, 5))) + + +def test_contrast_metric_bp2bpsim_rejects_invalid_bit_count() -> None: + with pytest.raises(ValueError, match="n_bits"): + compute_bp2bpsim(np.zeros((2, 2)), np.zeros((2, 2)), n_bits=0) diff --git a/tests/validation_tests/ImageEnhancement_validation_test.py b/tests/validation_tests/ImageEnhancement_validation_test.py new file mode 100644 index 0000000..48dcdbf --- /dev/null +++ b/tests/validation_tests/ImageEnhancement_validation_test.py @@ -0,0 +1,160 @@ +from pathlib import Path +import csv +import hashlib +import json + +import numpy as np +import pytest +from PIL import Image + +from shinier.color.Converter import rgb2gray +from shinier.utils import ( + MatlabOperators, + betce_gray, + compute_ambe, + compute_bp2bpsim, + compute_contrast_improvement, + compute_image_entropy, + compute_mssim, + compute_psnr, + nfldice_gray, + rdfhe_gray, + sfcef_gray, + tidhe_gray, +) + +pytestmark = pytest.mark.validation_tests + +# MATLAB reference hashes were generated from the free imTIDHE/imRDFHE/imNFLDICE/imBETCE/imSFCEF +# implementations linked in their respective articles; see tests/README.md. +ASSET_DIR = Path(__file__).resolve().parents[1] / "assets" +INPUT_DIR = ASSET_DIR / "SAMPLE_512X512" +SFCEF_MATLAB_DIR = ASSET_DIR / "sfcef_matlab_reference" +MATLAB_HASHES = json.loads((ASSET_DIR / "image_enhancement_matlab_sha256.json").read_text()) +ALGORITHMS = {"tidhe": tidhe_gray, "rdfhe": rdfhe_gray, "nfldice": nfldice_gray, "betce": betce_gray, "sfcef": sfcef_gray} +# SFCEF excluded: MATLAB filter2 uses FMA; pixel-level parity is tested separately via max_diff + exact diff count. +STRICT_HASH_ALGORITHMS = {"tidhe", "rdfhe", "nfldice", "betce"} +IMAGE_INDEX = {name: idx for idx, name in enumerate(MATLAB_HASHES["images"])} + + +def _matlab_gray(path: Path) -> np.ndarray: + rgb = np.asarray(Image.open(path).convert("RGB")) + return MatlabOperators.uint8(rgb2gray(rgb, weighting_standard="rec601", matlab_601=True)) + + +def _hash_output(algorithm: str, input_path: Path) -> dict[str, object]: + actual = ALGORITHMS[algorithm](_matlab_gray(input_path)) + expected_hash = MATLAB_HASHES["sha256"][algorithm][IMAGE_INDEX[input_path.name]] + return { + "algorithm": algorithm, + "image": input_path.name, + "shape": list(actual.shape), + "expected_shape": MATLAB_HASHES["shape"], + "dtype": str(actual.dtype), + "expected_dtype": MATLAB_HASHES["dtype"], + "sha256": hashlib.sha256(np.ascontiguousarray(actual).tobytes()).hexdigest(), + "expected_sha256": expected_hash, + } + + +@pytest.mark.parametrize("algorithm", STRICT_HASH_ALGORITHMS) +@pytest.mark.parametrize("input_path", sorted(INPUT_DIR.glob("*.png")), ids=lambda p: p.stem) +def test_image_enhancement_matches_matlab_reference_hash(algorithm: str, input_path: Path) -> None: + row = _hash_output(algorithm, input_path) + assert row["shape"] == row["expected_shape"], row + assert row["dtype"] == row["expected_dtype"], row + assert row["sha256"] == row["expected_sha256"], row + + +_SFCEF_SYNTHETIC_PATTERNS: dict[str, "np.ndarray"] = {} + + +def _sfcef_synthetic_inputs() -> dict[str, np.ndarray]: + """Three 1000×1000 synthetic patterns, all in [0, 32]. + + All pixels ≤ 32 guarantees 0% clipping in sfcef_gray (max output = 7.8×32 = 249.6 < 255), + exercising the full output range — unlike natural images where 80%+ of pixels saturate to 255. + Three gradient directions stress the filter in orthogonal orientations. + """ + if _SFCEF_SYNTHETIC_PATTERNS: + return _SFCEF_SYNTHETIC_PATTERNS + n = 1000 + ii, jj = np.mgrid[0:n, 0:n] + _SFCEF_SYNTHETIC_PATTERNS["sfcef_synth_diag"] = np.round(32.0 * (ii + jj) / (999 + 999)).astype(np.uint8) + _SFCEF_SYNTHETIC_PATTERNS["sfcef_synth_horiz"] = np.round(32.0 * jj / 999).astype(np.uint8) + _SFCEF_SYNTHETIC_PATTERNS["sfcef_synth_vert"] = np.round(32.0 * ii / 999).astype(np.uint8) + return _SFCEF_SYNTHETIC_PATTERNS + + +@pytest.mark.parametrize("pattern_name", ["sfcef_synth_diag", "sfcef_synth_horiz", "sfcef_synth_vert"]) +def test_sfcef_pixel_diff_from_matlab(pattern_name: str) -> None: + """SFCEF output on a synthetic low-value image must match MATLAB imSFCEF within FMA tolerance. + + Inputs: three 1000×1000 gradients in [0, 32] (zero clipping, full output range exercised). + Two assertions per pattern: + - max_diff <= 1: no pixel deviates by more than 1 gray level (FMA theoretical bound). + - n_diffs == expected: exact count of differing pixels, catching systematic shifts. + """ + gray = _sfcef_synthetic_inputs()[pattern_name] + matlab_ref = np.asarray(Image.open(SFCEF_MATLAB_DIR / f"{pattern_name}_matlab.png").convert("L")) + python_out = sfcef_gray(gray, legacy_mode=True) + diff = np.abs(python_out.astype(np.int16) - matlab_ref.astype(np.int16)) + n_diffs = int(np.sum(diff > 0)) + expected_n_diffs = MATLAB_HASHES["sfcef_synthetic_n_diffs"][pattern_name] + total_pixels = gray.size + print(f"\n {pattern_name}: {n_diffs}/{total_pixels} pixels off ({100 * n_diffs / total_pixels:.6f}%), max_diff={int(diff.max())}") + assert diff.max() <= 1, f"max pixel diff = {diff.max()} (expected ≤ 1)" + assert n_diffs == expected_n_diffs, ( + f"{n_diffs} differing pixels (expected {expected_n_diffs}). " + "A systematic shift would produce a much larger count." + ) + + +@pytest.mark.parametrize("pattern_name", ["sfcef_synth_diag", "sfcef_synth_horiz", "sfcef_synth_vert"]) +def test_sfcef_metric_diff_from_matlab_under_point_one_percent(pattern_name: str) -> None: + """SFCEF validation metrics must stay within 0.1% of MATLAB values. + + Pairwise metrics use the synthetic input as the reference image. CI and + entropy are computed directly on each enhanced output. + """ + gray = _sfcef_synthetic_inputs()[pattern_name] + matlab_ref = np.asarray(Image.open(SFCEF_MATLAB_DIR / f"{pattern_name}_matlab.png").convert("L")) + python_out = sfcef_gray(gray, legacy_mode=True) + metrics = { + "ambe_vs_input": (compute_ambe(gray, python_out), compute_ambe(gray, matlab_ref)), + "mssim_vs_input": (compute_mssim(gray, python_out, data_range=255), compute_mssim(gray, matlab_ref, data_range=255)), + "psnr_vs_input": (compute_psnr(gray, python_out, data_range=255), compute_psnr(gray, matlab_ref, data_range=255)), + "bp2bpsim_vs_input": (compute_bp2bpsim(gray, python_out), compute_bp2bpsim(gray, matlab_ref)), + "ci": (compute_contrast_improvement(python_out), compute_contrast_improvement(matlab_ref)), + "entropy": (compute_image_entropy(python_out), compute_image_entropy(matlab_ref)), + } + rows = {} + for name, (python_value, matlab_value) in metrics.items(): + abs_diff = float(np.abs(python_value - matlab_value)) + percent_of_matlab = 0.0 if matlab_value == 0 and abs_diff == 0 else abs_diff / float(np.abs(matlab_value)) * 100.0 + rows[name] = { + "python": float(python_value), + "matlab": float(matlab_value), + "abs_diff": abs_diff, + "percent_of_matlab": percent_of_matlab, + } + print(f"\n {pattern_name} metric diff vs MATLAB: {rows}") + assert all(row["percent_of_matlab"] <= 0.1 for row in rows.values()), rows + + +def test_image_enhancement_matlab_comparison_report(tmp_path: Path) -> None: + rows = [ + _hash_output(algorithm, path) + for algorithm in ALGORITHMS + for path in sorted(INPUT_DIR.glob("*.png")) + ] + report = tmp_path / "image_enhancement_matlab_python_hash_comparison.csv" + with report.open("w", newline="") as handle: + writer = csv.DictWriter(handle, fieldnames=list(rows[0])) + writer.writeheader() + writer.writerows(rows) + assert all( + row["sha256"] == row["expected_sha256"] + for row in rows + if row["algorithm"] in STRICT_HASH_ALGORITHMS + ), rows diff --git a/tests/validation_tests/ImageProcessor_validation_test.py b/tests/validation_tests/ImageProcessor_validation_test.py index 5fbe407..e74b548 100644 --- a/tests/validation_tests/ImageProcessor_validation_test.py +++ b/tests/validation_tests/ImageProcessor_validation_test.py @@ -1,4 +1,4 @@ -# tests/integration_tests/test_image_processor_validation_sharded.py +# tests/validation_tests/ImageProcessor_validation_test.py """Exhaustive ImageProcessor validations with smart pruning and sharding. Prunes only: @@ -13,30 +13,36 @@ Also restores RMSE–improvement checks for modes 5..8. Env: - SHARDS, SHARD_INDEX, SHOW_PROGRESS, DUMP_FILE_FORMAT, START_AT + COVERAGE_MODE, SHARDS, SHARD_INDEX, SHOW_PROGRESS, DUMP_FILE_FORMAT, + START_AT, PERCENT_SAMPLED (sampled mode only). See tests/README.md. """ from __future__ import annotations -import copy from traceback import format_exc -import itertools import os import shutil import random +import hashlib +import json from pathlib import Path from typing import Any, Callable, Dict, List, Optional, Tuple, Union, get_origin, get_args, Literal import numpy as np import pytest +from pydantic import ValidationError from tqdm.auto import tqdm -from shinier import ImageDataset, Options, utils, ImageListIO, ImageProcessor +from shinier import ImageDataset, Options, utils, ImageProcessor from tests import utils as utils_test from shinier.color.Converter import REC_STANDARD REC_STANDARD = [r for r in get_args(REC_STANDARD)] pytestmark = pytest.mark.validation_tests +# Exceptions Options() raises for *intentionally* invalid combos. Anything outside +# this set is treated as a real bug rather than silently classified as 'invalid'. +EXPECTED_INVALID = (ValidationError, ValueError, TypeError) + # -------------------- env -------------------- START_AT = int(os.getenv("START_AT", "0")) @@ -46,7 +52,12 @@ SHARDS = int(os.getenv("SHARDS", "1")) SHARD_INDEX = int(os.getenv("SHARD_INDEX", "0")) SHOW_PROGRESS = os.getenv("SHOW_PROGRESS", "1") == "1" -PERCENT_SAMPLED = float(os.getenv("PERCENT_SAMPLED", "1")) +PERCENT_SAMPLED = float(os.getenv("PERCENT_SAMPLED", "0.001")) +COVERAGE_MODE = os.getenv("COVERAGE_MODE", "sampled") # exhaustive | pruned | sampled +# Promote soft (quality-regression) anomalies to hard failures when set. +STRICT_SOFT = os.getenv("STRICT_SOFT", "0").lower() in ("1", "true", "yes") +# Aggregate gate: max fraction of executed combos allowed to regress (1.0 = off). +MAX_SOFT_FAIL_RATE = float(os.getenv("MAX_SOFT_FAIL_RATE", "1.0")) # SHARD_INDEX = 4 # SHARDS = 8 @@ -84,16 +95,13 @@ def test_imageprocessor_validations_sharded(test_tmpdir: Path) -> None: src_images_path = utils_test.get_small_imgs_path(utils_test.IMAGE_PATH) images_buffers = utils_test.prepare_images(utils_test.IMAGE_PATH) src0 = images_buffers['images'][0] - # with Image.open(src_images_path[0]) as pil_image: - # src0 = np.array(pil_image.convert("RGB")) h, w = src0.shape[:2] targets = utils_test.precompute_targets(images_buffers) - ag, ct, rs = 1, 0, 1 mask_dir = test_tmpdir / "MASK" utils_test.make_masks(mask_dir, h=h, w=w, n=1) - # ----- parameter grids (only pruned redundancies) ----- + # ----- parameter grids ----- choices = {name: get_possible_values(field) for name, field in Options.model_fields.items()} choices['input_folder'] = [utils_test.IMAGE_PATH] choices['masks_folder'] = [mask_dir] @@ -105,76 +113,141 @@ def test_imageprocessor_validations_sharded(test_tmpdir: Path) -> None: choices['hist_iterations'] = [3] choices['verbose'] = [-1] - all_fields = list(choices.keys()) - total_combo = np.prod([len(v) for v in choices.values() if hasattr(v, '__len__') and not isinstance(v, str)]) # Some of which are not valid and duplicated + # Pruned grid: collapses options whose values share the same code path. + # Collapsed values are individually verified in PrunedCoverage_test.py. + choices_pruned = dict(choices) + choices_pruned['fft_padding_mode'] = [0, 1] # off vs. on; modes 2/3 are non-distinct + choices_pruned['standalone_op'] = ['ie_methods'] # dithering path covered by dithering param above + choices_pruned['ie_methods'] = ['tidhe'] # all methods share the same call site + choices_pruned['rec_standard'] = [1] # same code path, different matrices only + choices_pruned['seed'] = [None] # only affects tie-breaking randomness + choices_pruned['dithering'] = [0, 1] # off + one ordered method + choices_pruned['hist_specification'] = [1, 3] # basic + one advanced variant + choices_pruned['rescaling'] = [0, 2] # off + default; variants 1/3 non-distinct + choices_pruned['conserve_memory'] = [False] # no code-path interaction with other params + choices_pruned['background'] = [300, 120] # 300 = automatic (special sentinel), 120 = manual + + active_choices = choices_pruned if COVERAGE_MODE == "pruned" else choices + all_fields = list(active_choices.keys()) + sizes = [len(active_choices[f]) for f in all_fields] + total_combo = int(np.prod(sizes)) + + # Fingerprint of the option space (excluding volatile per-run paths). When the + # grid changes, the fingerprint changes, so cached index->combo entries from an + # incompatible space are simply treated as not-done instead of mis-attributed. + _volatile = ("input_folder", "masks_folder", "output_folder") + schema_fp = hashlib.sha1( + json.dumps( + [(f, [str(v) for v in active_choices[f]]) for f in all_fields if f not in _volatile], + sort_keys=False, + ).encode() + ).hexdigest()[:12] # Initialize the db if RESTART: - utils_test.initialize_db(total_number_of_tests=total_combo) + utils_test.initialize_db() if START_AT: - utils_test.mark_hash_range_done(start=0, end=START_AT) + utils_test.mark_hash_range_done(start=0, end=START_AT, namespace=COVERAGE_MODE) + + def _index_to_combo(flat_idx: int) -> tuple: + """Decode a flat combo index into per-field indices (row-major).""" + indices = [] + for s in reversed(sizes): + indices.append(flat_idx % s) + flat_idx //= s + return tuple(reversed(indices)) + + # Generate indices to visit. + # Sampled: draw random indices directly (O(n_to_draw)), avoiding O(total_combo) traversal. + # Exhaustive/pruned: stride so each shard processes exactly 1/SHARDS of all combos. + if COVERAGE_MODE == "sampled": + # Set PRNG seed per shard for reproducible, independent samples + rng = random.Random(int(f'98234987234{SHARD_INDEX}')) + n_to_draw = max(1, int(total_combo * PERCENT_SAMPLED / max(SHARDS, 1))) + if n_to_draw < total_combo // 2: + indices_to_process: Any = sorted(rng.sample(range(total_combo), n_to_draw)) + else: + indices_to_process = sorted(set(rng.randint(0, total_combo - 1) for _ in range(n_to_draw))) + n_indices = len(indices_to_process) + else: + indices_to_process = range(SHARD_INDEX, total_combo, SHARDS) + n_indices = len(indices_to_process) pbar = None if SHOW_PROGRESS and tqdm is not None: - per_shard = total_combo // SHARDS + (1 if SHARD_INDEX < (total_combo % SHARDS) else 0) - pbar = tqdm(total=total_combo, initial=START_AT, desc=f"Shard {SHARD_INDEX+1}/{SHARDS}", ncols=0) - - # Set PNRG seed for test sampling - rng = random.Random() # independent state - rng.seed(int(f'98234987234{SHARD_INDEX}')) - - # Reset the iterator for the real loop + pbar = tqdm(total=n_indices, initial=0, desc=f"[{COVERAGE_MODE}] Shard {SHARD_INDEX+1}/{SHARDS}", ncols=0) + + failures: list[str] = [] # hard failures: unexpected exceptions → test fails + soft_failures: list[str] = [] # soft failures: reported but test passes + # Coverage accounting (#7): every combo past START_AT lands in exactly one bucket. + executed = 0 # ran the pipeline successfully + skipped_cached = 0 # skipped because already terminal-success in the DB + invalid = 0 # intentionally invalid option combination + hard_failed = 0 # exception / dumped hard failure + regressed = 0 # executed but produced soft anomalies cnt = -1 - for i, combo in enumerate(itertools.product(*(choices[f] for f in all_fields))): + for i in indices_to_process: if pbar is not None: pbar.update(1) - combo_hash = str(i) - kwargs = dict(zip(all_fields, combo)) - - # if `i` within this SHARD_INDEX, proceed else next `i` - if i % SHARDS != SHARD_INDEX: + # if `i` >= START_AT, proceed else next `i` (exhaustive/pruned only) + if COVERAGE_MODE != "sampled" and i < START_AT: continue - cnt += 1 - # if `i` >= START_ITER, proceed else next `i` - if i < START_AT: - continue + # Namespaced by coverage-mode + option-space fingerprint so neither the + # different modes nor an evolving grid collide in the shared DB. + combo_hash = f"{COVERAGE_MODE}:{schema_fp}:{i}" + combo = tuple(active_choices[f][idx] for f, idx in zip(all_fields, _index_to_combo(i))) + kwargs = dict(zip(all_fields, combo)) + cnt += 1 # if combo never tested, proceed else next `i` if utils_test.is_already_done(combo_hash): - continue - - # if not randomly selected, proceed else next `i` - is_sampled = rng.random() >= (1 - PERCENT_SAMPLED) - if not is_sampled: - # Register that combo as invalid - utils_test.mark_hash_status(combo_hash, status='not_sampled') + skipped_cached += 1 continue # Create target hist or spectrum if requested if (kwargs['mode'] == 9 and kwargs['dithering'] == 0) or (kwargs['mode'] == 9 and kwargs['legacy_mode']): # Register that combo as invalid + invalid += 1 utils_test.mark_hash_status(combo_hash, status='invalid_option_combination') continue try: opts = Options(**kwargs) + except EXPECTED_INVALID as e: + # Expected validation rejection → legitimately invalid combo + invalid += 1 + utils_test.mark_hash_status(combo_hash, status='invalid', error=get_error_msg(e)) + continue except Exception as e: - # Register that combo as invalid and set error as message from exception - error_msg = get_error_msg(e) - utils_test.mark_hash_status(combo_hash, status='invalid', error=error_msg) + # Unexpected exception type while building Options → real bug, not "invalid" + hard_failed += 1 + failures.append( + f"\n💥 Options() raised unexpected {e.__class__.__name__} (combo {i})\n" + f"→ Combo: {kwargs}\n" + f"→ {e}\n" + ) + utils_test.mark_hash_status(combo_hash, status='error', error=get_error_msg(e)) continue ag = True if kwargs['legacy_mode'] else kwargs['as_gray'] opts.target_hist = targets["hist"][int(ag)][kwargs['linear_luminance']][kwargs['rec_standard']] if kwargs['target_hist'] == "unit_test" else ("equal" if kwargs['target_hist'] == "equal" else None) - opts.target_spectrum = targets["spec"][int(ag)][kwargs['linear_luminance']][kwargs['rec_standard']] if kwargs['target_spectrum'] == "unit_test" else None + # Precomputed spectra are at original image size; fft_padding_mode changes the expected + # size, so let the pipeline compute its own target when padding is active. + if kwargs['target_spectrum'] == 'unit_test' and kwargs.get('fft_padding_mode', 0) != 0: + opts.target_spectrum = None + else: + opts.target_spectrum = targets["spec"][int(ag)][kwargs['linear_luminance']][kwargs['rec_standard']] if kwargs['target_spectrum'] == "unit_test" else None - # Set seed + # Set seed. Derive it from a combo with volatile per-run path fields stripped + # (e.g. the tmp masks_folder), otherwise the seed — and thus the pipeline output — + # would change every run and reintroduce flakiness. if opts.seed is not None: seed_iter = opts.seed else: - seed_iter = utils_test.deterministic_seed_from_combo(combo=combo) + stable_combo = tuple(v for f, v in zip(all_fields, combo) if f not in _volatile) + seed_iter = utils_test.deterministic_seed_from_combo(combo=stable_combo) # Log file out_dir = test_tmpdir / ( @@ -193,94 +266,140 @@ def test_imageprocessor_validations_sharded(test_tmpdir: Path) -> None: f"_lm{int(opts.legacy_mode)}" # Legacy mode ) out_dir.mkdir(parents=True, exist_ok=True) - ag2, ct, rs = None, None, None try: opts.output_folder = out_dir - ag2, ct, rs = int(opts.as_gray), int(opts.linear_luminance), opts.rec_standard except Exception as e: + invalid += 1 if out_dir and out_dir.exists(): shutil.rmtree(out_dir, ignore_errors=True) # Register that combo as invalid and set error as message from exception - error_msg = get_error_msg(e) - utils_test.mark_hash_status(combo_hash, status='invalid', error=error_msg) + utils_test.mark_hash_status(combo_hash, status='invalid', error=get_error_msg(e)) continue rand_selected_paths = None + combo_warnings: list[str] = [] try: # Options / pipeline rand_selected_images = utils_test.select_n_imgs(images_buffers['images'], n=2, seed=seed_iter) # sRGB - rand_selected_buffers = utils_test.select_n_imgs(images_buffers['buffers'][ag2][ct][rs], n=2, seed=seed_iter) # Y from xyY or sRGB (depending on linear_luminance) rand_selected_paths = utils_test.select_n_imgs(src_images_path, n=2, seed=seed_iter) # sRGB - # as_gray_ds = 1 if opts.as_gray == True and opts.linear_luminance is True else 0 - - # images_copy = ImageListIO(input_data=rand_selected_images, conserve_memory=False) - initial_buffers = ImageListIO(input_data=rand_selected_buffers, conserve_memory=False) - # buffers_empty = [np.zeros(im.shape, dtype=bool) for im in initial_buffers] ds = ImageDataset.model_construct(images=rand_selected_images, options=opts) - proc = ImageProcessor.model_construct(dataset=ds, options=opts, verbose=-1, from_validation_test=True) + # from_unit_test=True makes post_init skip the auto-run so the seed can be set + # deterministically before processing; from_validation_test=True keeps the SSIM + # validation recording enabled. This makes every combo reproducible and ensures + # the dumped seed is the one actually used. + proc = ImageProcessor.model_construct( + dataset=ds, options=opts, verbose=-1, + from_unit_test=True, from_validation_test=True, + ) + proc.seed = seed_iter + proc.process() # Prepare targets for validation th = opts.target_hist if proc._target_hist is None else proc._target_hist ts = opts.target_spectrum if proc._target_spectrum is None else proc._target_spectrum - # Prepare images for validation: convert them into xyY if needed + # Use the processor's own buffers so the before/after comparison shares the + # exact color domain the pipeline produced (forward color treatment included). + initial_buffers = proc._initial_buffer final_buffers = proc._final_buffer - # internal validations + # internal validations — hard fail (crash = real bug) for rec in getattr(proc, "validation", []): - if not bool(rec.get("valid_result", True)): + if _is_fail(rec.get("valid_result", True)): _dump_and_fail(rec, kwargs, seed_iter, rand_selected_paths, test_tmpdir) - # --------- RESTORED: RMSE improvement checks for modes 5..8 ---------- + # --------- RMSE improvement checks for modes 5..8 ---------- + # Soft fail (strict): RMSE must decrease (improve); any regression is noted. + # Hard fail (lenient): RMSE more than doubled or increased > 0.1 — something badly wrong. + _rmse_hard_tol = lambda b: np.maximum(0.1, b * 2) + 1e-9 # noqa: E731 if opts.mode in (5, 6, 7, 8): # Histogram (always for 5..8) _, rmse_hist_before = utils.hist_match_validation(images=initial_buffers, binary_masks=proc.bool_masks, target_hist=th, normalize_rmse=True) _, rmse_hist_after = utils.hist_match_validation(images=final_buffers, binary_masks=proc.bool_masks, target_hist=th, normalize_rmse=True) - if not np.all(rmse_hist_after + 1e-9 <= rmse_hist_before): - rec = { - "iter": -1, "step": -1, "processing_function": "hist_match", - "valid_result": False, - "log_result": f"Histogram RMSE not improved: {rmse_hist_before} -> {rmse_hist_after}", - } - _dump_and_fail(rec, kwargs, seed_iter, rand_selected_paths, test_tmpdir) + if not np.all(rmse_hist_after <= rmse_hist_before + _rmse_hard_tol(rmse_hist_before)): + _dump_and_fail({"iter": -1, "step": -1, "processing_function": "hist_match", "valid_result": False, + "log_result": f"Histogram RMSE more than doubled: {rmse_hist_before} -> {rmse_hist_after}"}, + kwargs, seed_iter, rand_selected_paths, test_tmpdir) + elif not np.all(rmse_hist_after + 1e-9 <= rmse_hist_before): + combo_warnings.append(f"Histogram RMSE not improved: {rmse_hist_before} -> {rmse_hist_after}") # Spatial frequency (5,7) if opts.mode in (5, 7): - _, rmse_sf_before = utils.sf_match_validation(images=initial_buffers, target_spectrum=ts, normalize_rmse=True) - _, rmse_sf_after = utils.sf_match_validation(images=final_buffers, target_spectrum=ts, normalize_rmse=True) - if not np.all(rmse_sf_after <= rmse_sf_before + 1e-9): - rec = { - "iter": -1, "step": -1, "processing_function": "sf_match", - "valid_result": False, - "log_result": f"SF RMSE not improved: {rmse_sf_before} -> {rmse_sf_after}", - } - _dump_and_fail(rec, kwargs, seed_iter, rand_selected_paths, test_tmpdir) + _, rmse_sf_before = utils.sf_match_validation(images=initial_buffers, target_spectrum=ts, normalize_rmse=True, fft_padding_mode=opts.fft_padding_mode, fft_padding_value=opts.fft_padding_value) + _, rmse_sf_after = utils.sf_match_validation(images=final_buffers, target_spectrum=ts, normalize_rmse=True, fft_padding_mode=opts.fft_padding_mode, fft_padding_value=opts.fft_padding_value) + if not np.all(rmse_sf_after <= rmse_sf_before + _rmse_hard_tol(rmse_sf_before)): + _dump_and_fail({"iter": -1, "step": -1, "processing_function": "sf_match", "valid_result": False, + "log_result": f"SF RMSE more than doubled: {rmse_sf_before} -> {rmse_sf_after}"}, + kwargs, seed_iter, rand_selected_paths, test_tmpdir) + elif not np.all(rmse_sf_after + 1e-9 <= rmse_sf_before): + combo_warnings.append(f"SF RMSE not improved: {rmse_sf_before} -> {rmse_sf_after}") # Spectrum (6,8) if opts.mode in (6, 8): - _, rmse_spec_before = utils.spec_match_validation(images=initial_buffers, target_spectrum=ts, normalize_rmse=True) - _, rmse_spec_after = utils.spec_match_validation(images=final_buffers, target_spectrum=ts, normalize_rmse=True) - if not np.all(rmse_spec_after <= rmse_spec_before + 1e-9): - rec = { - "iter": -1, "step": -1, "processing_function": "spec_match", - "valid_result": False, - "log_result": f"Spectrum RMSE not improved: {rmse_spec_before} -> {rmse_spec_after}", - } - _dump_and_fail(rec, kwargs, seed_iter, rand_selected_paths, test_tmpdir) + _, rmse_spec_before = utils.spec_match_validation(images=initial_buffers, target_spectrum=ts, normalize_rmse=True, fft_padding_mode=opts.fft_padding_mode, fft_padding_value=opts.fft_padding_value) + _, rmse_spec_after = utils.spec_match_validation(images=final_buffers, target_spectrum=ts, normalize_rmse=True, fft_padding_mode=opts.fft_padding_mode, fft_padding_value=opts.fft_padding_value) + if not np.all(rmse_spec_after <= rmse_spec_before + _rmse_hard_tol(rmse_spec_before)): + _dump_and_fail({"iter": -1, "step": -1, "processing_function": "spec_match", "valid_result": False, + "log_result": f"Spectrum RMSE more than doubled: {rmse_spec_before} -> {rmse_spec_after}"}, + kwargs, seed_iter, rand_selected_paths, test_tmpdir) + elif not np.all(rmse_spec_after + 1e-9 <= rmse_spec_before): + combo_warnings.append(f"Spectrum RMSE not improved: {rmse_spec_before} -> {rmse_spec_after}") + + # Recognizability safeguard on the FINAL composite output (not per-step + # monotonicity). Two thresholds, mirroring the RMSE gate above: + # - hard (collapse): the image is structurally destroyed -> fail. + # - soft (degraded): still recognizable but quality is dropping -> warn. + # Composite SSIM legitimately falls as iterations rise (the hist<->spectrum + # constraints are incompatible and fight each other), so the hard floor stays + # low to catch only a destroyed image, while the soft floor surfaces the + # degraded regime where real problems begin. Both are env-tunable. + ssim_collapse = float(os.getenv("COMPOSITE_SSIM_FLOOR", "0.3")) + ssim_warn = float(os.getenv("COMPOSITE_SSIM_WARN", "0.6")) + ssim_means = [] + for _idx in range(len(final_buffers)): + _, _ssim = utils.ssim_sens(initial_buffers[_idx], final_buffers[_idx], + data_range=255, use_sample_covariance=False, + binary_mask=proc.bool_masks[_idx]) + ssim_means.append(float(np.mean(_ssim))) + mean_ssim = float(np.mean(ssim_means)) + if mean_ssim < ssim_collapse: + _dump_and_fail({"iter": -1, "step": -1, "processing_function": "composite", "valid_result": False, + "log_result": f"Final composite SSIM {mean_ssim:.4f} collapsed below {ssim_collapse} (per-image {np.round(ssim_means, 4)})"}, + kwargs, seed_iter, rand_selected_paths, test_tmpdir) + elif mean_ssim < ssim_warn: + combo_warnings.append(f"Composite SSIM degraded: {mean_ssim:.4f} < {ssim_warn} (per-image {np.round(ssim_means, 4)})") # --------------------------------------------------------------------- - # SSIM optimization monotonicity + # SSIM optimization checks — per-hist_match passes. + # The 'final' pass checks SSIM(first proposal -> hist_match output). This is a + # meaningful *terminal* check only for mode 2, where hist_match is the last op. + # In composite modes (5-8) hist_match is intermediate — sf_match/spec_match + # deliberately move the image afterwards — so per-step SSIM monotonicity is not + # psychophysically meaningful; it is demoted to diagnostic and the composite + # end-state SSIM floor below is the real recognizability gate. if opts.hist_optim and opts.mode in (2, 5, 6, 7, 8): - ssim_records = getattr(proc, "ssim_results", []) - if ssim_records and not all(bool(r.get("valid_result", True)) for r in ssim_records): - rec = {"iter": -1, "step": -1, "processing_function": "hist_match", - "valid_result": False, "log_result": f"SSIM optimization non-monotonic: {ssim_records}"} - _dump_and_fail(rec, kwargs, seed_iter, rand_selected_paths, test_tmpdir) - except ControlledFailure: - pass # handled already + for rec in getattr(proc, "ssim_results", []): + if _is_fail(rec.get("valid_result", True)): + tag = rec.get('tag', 'sub_iter') + msg = ( + f"SSIM optimization {tag} regression: " + f"image={rec.get('image')}, channel={rec.get('channel')}, " + f"iter={rec.get('iter')}, step={rec.get('step')}, " + f"ssim_values={rec.get('ssim_values')}" + ) + if tag == 'final' and opts.mode == 2: + _dump_and_fail({**rec, "log_result": msg}, kwargs, seed_iter, rand_selected_paths, test_tmpdir) + else: + combo_warnings.append(msg) + + except ControlledFailure as cf: + hard_failed += 1 + failures.append(str(cf)) + utils_test.mark_hash_status(combo_hash, status='failed') except Exception as e: + hard_failed += 1 tb = format_exc() dump_path = utils_test.dump_failure_context( combo_dict=kwargs, @@ -290,26 +409,96 @@ def test_imageprocessor_validations_sharded(test_tmpdir: Path) -> None: selected_paths=rand_selected_paths, file_type=DUMP_FILE_FORMAT, ) - raise AssertionError( + failures.append( f"\n💥 Unexpected error\n" f"→ Shard {SHARD_INDEX}, Combo global-index {i}, Per-shard #{cnt}\n" f"→ Combo: {kwargs}\n" f"→ Exception: {e.__class__.__name__}: {e}\n" f"→ Dumped context: {dump_path}\n" ) + utils_test.mark_hash_status(combo_hash, status='error') else: + # Pipeline ran to completion for this combo. + executed += 1 # Cleanup safely if out_dir and out_dir.exists(): shutil.rmtree(out_dir, ignore_errors=True) - # Register that combo as done - utils_test.mark_hash_status(combo_hash, status='done') + if combo_warnings: + regressed += 1 + warn_str = "; ".join(combo_warnings) + dump_path = utils_test.dump_failure_context( + combo_dict=kwargs, + rec={"iter": -1, "step": -1, "warnings": combo_warnings}, + tmp_root=test_tmpdir, + seed=seed_iter, + selected_paths=rand_selected_paths or [], + file_type=DUMP_FILE_FORMAT, + ) + entry = ( + f"\nSoft failure — Combo {i} mode={opts.mode} hist_optim={opts.hist_optim}\n" + f" {warn_str}\n" + f" → Dumped: {dump_path}" + ) + if STRICT_SOFT: + # Promote to a hard failure and keep it red on re-run until fixed. + failures.append(entry) + utils_test.mark_hash_status(combo_hash, status='failed', error=warn_str) + else: + soft_failures.append(entry) + utils_test.mark_hash_status(combo_hash, status='done', error=warn_str) + else: + utils_test.mark_hash_status(combo_hash, status='done') if pbar is not None: pbar.close() - assert cnt >= 0 + if soft_failures: + print(f"\n[soft failures] {len(soft_failures)} combo(s) reported optimization anomalies:") + for sf in soft_failures: + print(f" {sf}") + + # ----- coverage accounting (#7) ----- + # Every combo that passed the START_AT filter (cnt + 1 of them) must land in + # exactly one bucket; otherwise a combo was silently dropped. + total_seen = cnt + 1 + bucketed = executed + skipped_cached + invalid + hard_failed + print( + f"\n[coverage] seen={total_seen} executed={executed} cached={skipped_cached} " + f"invalid={invalid} hard_failed={hard_failed} regressed={regressed}" + ) + assert bucketed == total_seen, ( + f"combo accounting mismatch: {bucketed} bucketed != {total_seen} seen " + "(a combo was silently dropped)" + ) + # Prove the shard did real work this run (ran/evaluated at least one combo) — + # unless every combo was legitimately cached from a previous run. + assert (executed + hard_failed + invalid) > 0 or skipped_cached == total_seen, ( + "no combos were evaluated and not all were cached — the run did no real work" + ) + + # ----- aggregate soft-failure gate (#5) ----- + soft_rate = regressed / max(executed, 1) + assert soft_rate <= MAX_SOFT_FAIL_RATE, ( + f"soft-failure (regression) rate {soft_rate:.2%} exceeds " + f"MAX_SOFT_FAIL_RATE={MAX_SOFT_FAIL_RATE:.2%}" + ) + + assert not failures, f"{len(failures)} combo(s) failed:\n" + "\n---\n".join(failures) # -------------------- helpers -------------------- +def _is_fail(valid_result: Any) -> bool: + """Normalize the mixed ``valid_result`` field into a failure boolean. + + ``valid_result`` can be a bool / numpy bool (from ``_validate``) or one of the + strings ``"PASS"`` / ``"WARN"`` / ``"FAIL"`` (from diagnostic records). Only + a boolean-False or the explicit ``"FAIL"`` string counts as a failure; + ``"WARN"`` is a diagnostic and must not fail the test. + """ + if isinstance(valid_result, str): + return valid_result == "FAIL" + return not bool(valid_result) + + class ControlledFailure(AssertionError): """Raised when _dump_and_fail already handled the failure context."""