diff --git a/docs/apidoc.sh b/docs/apidoc.sh index a25c6152..f8b7fb19 100644 --- a/docs/apidoc.sh +++ b/docs/apidoc.sh @@ -27,7 +27,7 @@ echo "Updating generated RST files..." # sed -i '0,/dataset.dataset/s/dataset.dataset/pyvale.dataset/' source/pyvale.dataset.dataset.rst || error_exit "Failed to update dataset title" # Modules to process -modules=("dic" "blender" "mooseherder" "sensorsim" "verif" "dataset") +modules=("dic" "blender" "mooseherder" "sensorsim" "verif" "dataset" "specklegen") for mod in "${modules[@]}"; do rst="source/pyvale.${mod}.rst" diff --git a/docs/source/api_py.rst b/docs/source/api_py.rst index 1600d894..1f3f01e0 100644 --- a/docs/source/api_py.rst +++ b/docs/source/api_py.rst @@ -13,4 +13,6 @@ Detailed Python API pyvale.mooseherder pyvale.verif pyvale.dataset + pyvale.specklegen + diff --git a/docs/source/conf.py b/docs/source/conf.py index ffe15bee..c8e5ea77 100644 --- a/docs/source/conf.py +++ b/docs/source/conf.py @@ -117,6 +117,7 @@ '../../src/pyvale/examples/dic', '../../src/pyvale/examples/blenderimagedef', '../../src/pyvale/examples/mooseherder', + '../../src/pyvale/examples/specklegen', ], # Path to where to save gallery generated output 'gallery_dirs': [ @@ -125,6 +126,7 @@ 'examples/dic', 'examples/blenderimagedef', 'examples/mooseherder', + 'examples/specklegen', ], # Pattern to identify example files 'filename_pattern': '/plot_', diff --git a/docs/source/examples/examples.rst b/docs/source/examples/examples.rst index f3f38c02..2f077375 100644 --- a/docs/source/examples/examples.rst +++ b/docs/source/examples/examples.rst @@ -37,3 +37,10 @@ Mooseherder :maxdepth: 2 examples_mooseherder + +Specklegen +------------------------------------- +.. toctree:: + :maxdepth: 2 + + examples_specklegen diff --git a/docs/source/examples/examples_specklegen.rst b/docs/source/examples/examples_specklegen.rst new file mode 100644 index 00000000..5b22f9ec --- /dev/null +++ b/docs/source/examples/examples_specklegen.rst @@ -0,0 +1,16 @@ +.. _examples_specklegen: + +Specklegen +================= + +.. toctree:: + :maxdepth: 1 + + specklegen/ex1a_random_disks_overlap.rst + specklegen/ex1b_random_disks_reduce_overlap.rst + specklegen/ex1c_random_disks_grid.rst + specklegen/ex2a_perlin_noise.rst + specklegen/ex3a_fractal_noise.rst + specklegen/ex4a_simplex_noise.rst + + diff --git a/pyproject.toml b/pyproject.toml index d3947d81..d642a692 100644 --- a/pyproject.toml +++ b/pyproject.toml @@ -33,9 +33,13 @@ dependencies = [ "pybind11>=2.13.6", "pyqtgraph>=0.13.7", "opencv-python<=4.9.0.80", + "seaborn>=0.13.2", + "scikit-image>=0.0", + "perlin_numpy>=0.0.0", + "opensimplex>=0.4.5.1", "pandas>=2.3.1", "scikit-build-core>=0.11.6", - "ninja>=1.13.0" + "ninja>=1.13.0", ] [tool.pytest.ini_options] diff --git a/scripts/gengold_specklegen.py b/scripts/gengold_specklegen.py new file mode 100644 index 00000000..97504b37 --- /dev/null +++ b/scripts/gengold_specklegen.py @@ -0,0 +1,42 @@ +#=============================================================================== +# pyvale: the python validation engine +# License: MIT +# Copyright (C) 2025 The Computer Aided Validation Team +#=============================================================================== + +import pyvale.verif.specklegold as specklegold +import pyvale.verif.specklegenconst as specklegenconst + +def main() -> None: + + tags = ["random_disks", "random_disks_grid", "perlin", "fractal", "simplex"] + + for tag in tags: + + param_dict = { + "speckle_size": 20, + "screen_size_width": 1000, + "screen_size_height": 800, + "bit_depth": 8, + "theme": 'white_on_black', + "seed": 123, + "type_gen": tag, + "octaves": 3, + "lacunarity": 2, + "sigma": 4.0, + "reduce_overlap": True, + "attempts_tot": 300, + "perturbation_max": 12, + "case_tot": 3 + } + + print(80*"=") + print(f"Gold Output Generator for pyvale {tag} speckle pattern generation") + print(80*"=") + print(f"Saving gold output to: {specklegenconst.GOLD_PATH}\n") + + print(f"Generating gold output for {tag} field point sensors...") + specklegold.gen_gold_measurements(param_dict) + +if __name__ == "__main__": + main() diff --git a/src/pyvale/__init__.py b/src/pyvale/__init__.py index 371a5c8c..430adadd 100644 --- a/src/pyvale/__init__.py +++ b/src/pyvale/__init__.py @@ -21,3 +21,4 @@ from . import mooseherder from . import dataset from . import calib +from . import specklegen diff --git a/src/pyvale/examples/specklegen/README.md b/src/pyvale/examples/specklegen/README.md new file mode 100644 index 00000000..fbbd16fc --- /dev/null +++ b/src/pyvale/examples/specklegen/README.md @@ -0,0 +1,2 @@ +Specklegen Examples +================== diff --git a/src/pyvale/examples/specklegen/ex1a_random_disks_overlap.py b/src/pyvale/examples/specklegen/ex1a_random_disks_overlap.py new file mode 100644 index 00000000..84e3a265 --- /dev/null +++ b/src/pyvale/examples/specklegen/ex1a_random_disks_overlap.py @@ -0,0 +1,142 @@ +# ============================================================================== +# pyvale: the python validation engine +# License: MIT +# Copyright (C) 2025 The Computer Aided Validation Team +# ============================================================================== + +""" +Specklegen: Speckle pattern generation using random disk placement without checking for overlap +================================================================================ +Script to generate a synthetic speckle pattern made from randomly placed circular +speckles (disks), run diagnostics on the generated image, and save both the image +and diagnostics to the selected folder. +""" + +import numpy as np +import time +import json +import os +import pyvale.specklegen as specklegen + +#%% +# Here we parse command line arguments to set the speckle pattern parameters. +# For ease of use in this example script we set parameter values directly in the +# code rather than via bash script. +# The parameter responsible for reducing overlap is set to 'False' in this example. + +speckle_size = 20 +screen_size_width = 1000 +screen_size_height = 800 +bit_depth = 8 +theme = 'white_on_black' +seed = 10 +sigma = 4.0 +reduce_overlap = False +type_gen = "random_disks" +output_path = "src/pyvale/examples/specklegen/output/ex1a" + +print('Start') + +assert theme in ['black_on_white', 'white_on_black'], "Theme should be either 'black_on_white' or 'white_on_black'." + +if reduce_overlap: + print("Reducing overlap between speckles") +else: + print("Not reducing overlap between speckles") + +subfolder = f"/{type_gen}_{speckle_size}_{screen_size_width}_{screen_size_height}_{bit_depth}_{theme}_{seed}_{sigma}_{reduce_overlap}" +print(subfolder) +save_path = output_path + subfolder +if not os.path.exists(save_path): + os.makedirs(save_path) + +#%% +# We calculate parameteres aiming for approximately 50/50 black-to-white ratio. +# We now generate the speckle pattern using the specified parameters. +# The background and foreground colours are set based on the chosen theme and bit depth. + +speckle_area = np.pi * (speckle_size / 2) ** 2 +total_area = screen_size_width * screen_size_height +total_speckles = int((0.5 * total_area) / speckle_area) +print(f"Total number of speckles = {total_speckles}") +dynamic_range: int = 2**bit_depth - 1 +background_colour = 0 if theme == 'white_on_black' else dynamic_range +foreground_colour = dynamic_range if theme == 'white_on_black' else 0 + +feature_size_width = speckle_size +feature_size_height = speckle_size + +time_start = time.time() +image, results = specklegen.generate_speckles(screen_size_width, screen_size_height, + feature_size_width, feature_size_height, + foreground_colour, background_colour, + bit_depth, type_gen, seed, + total_speckles=total_speckles, + reduce_overlap=reduce_overlap, + sigma=sigma) +time_end = time.time() +time_taken = time_end - time_start +print(f"Time taken for speckle generation: {np.round(time_taken, 3)} seconds") + +np.savetxt(f"{save_path}/speckle_placement_results.csv", results, delimiter=",", + header="speckle_number, attempts, overlap(1/0/2), cent_x, cent_y", comments='', fmt=['%d', '%d', '%d', '%.3f', '%.3f']) + +#%% +# Now we run diagnostics on the generated speckle pattern and save the results. +# Finally, we print out the key statistics to the console. +# The plots are saved in the provided output folder. However, the diagnostic function outputs the matplotlib figures and axes, +# so the plot formatting could be changed from the default one used by the function. +# We aim to achieve black-to-white ratio as close to unity as possible. Unity ratio means 50/50 distribution of black and white colours. +# However, in this example, black-to-white ratio considerably deviates from unity. +# It is also visible in the irradiance value histogram. +# The proportion of 0 irradiance values, corresponding to black colour, overweighs the 255 values, corresponding to white colour. +# This is a result of speckle overlap. + +print("") +print('Starting speckle pattern diagnostics...') +results = specklegen.speckle_pattern_statistics(image, dynamic_range) +plots = specklegen.speckle_pattern_plots(image, dynamic_range, save_path) + +with open(f"{save_path}/speckle_pattern_diagnostics.json", 'w') as f: + json.dump(results, f, indent=4) + +ratio = results.get("black_white_ratio", None) +mean_gradient = results.get("mean_intensity_gradient", None) +std_dev = results.get("std_dev_irradiance", None) +avg = results.get("avg_irradiance", None) +contrast = results.get("contrast", None) +entropy = results.get("shannon_entropy", None) +peak_to_mean = results.get("peak_to_mean_ratio", None) +skew = results.get("skewness", None) +kurt = results.get("kurtosis", None) +avg_speckle_size_fwhm = results.get("avg_speckle_size_fwhm", None) +avg_speckle_size_e2 = results.get("avg_speckle_size_e2", None) +H_fit_stats = results.get("H_fit_stats", None) +V_fit_stats = results.get("V_fit_stats", None) + +print("") +print("Speckle statistics:") + +print(f"Black/White ratio: {np.round(ratio, 3)}") +print(f"Mean intensity gradient: {np.round(mean_gradient, 3)}") +print(f"Standard deviation of irradiance values: {np.round(std_dev, 3)}") +print(f"Average irradiance value: {np.round(avg, 3)}") +print(f"Contrast (std/mean): {np.round(contrast, 3)}") +print(f"Skewness: {np.round(skew, 3)}") +print(f"Kurtosis: {np.round(kurt, 3)}") +print(f"Shannon entropy: {np.round(entropy, 3)}") +print(f"Peak to mean ratio: {np.round(peak_to_mean, 3)}") +print(f"Average speckle size (full width at half maximum): {np.round(avg_speckle_size_fwhm, 3)} pixels") +print(f"Average speckle size (1/e^2): {np.round(avg_speckle_size_e2, 3)} pixels") +print(f"R_squared: Horisontal fit: {np.round(H_fit_stats['R_squared'], 3)}, Vertical fit: {np.round(V_fit_stats['R_squared'], 3)}") + +#%% +# Finally, the relative errors beetween the specified speckle size and the speckle size approximated using cautocovariance are calculated. +error = np.abs(avg_speckle_size_fwhm - speckle_size) * 100 / speckle_size +print(f"Percentage error between requested speckle size and measured speckle size from FWHM: {np.round(error, 3)} %") +error = np.abs(avg_speckle_size_e2 - speckle_size) * 100 / speckle_size +print(f"Percentage error between requested speckle size and measured speckle size from 1/e^2: {np.round(error, 3)} %") +np.save(f"{save_path}/image.npy", image) +print("") +print('End :)') +print("") \ No newline at end of file diff --git a/src/pyvale/examples/specklegen/ex1b_random_disks_reduce_overlap.py b/src/pyvale/examples/specklegen/ex1b_random_disks_reduce_overlap.py new file mode 100644 index 00000000..d5f31dea --- /dev/null +++ b/src/pyvale/examples/specklegen/ex1b_random_disks_reduce_overlap.py @@ -0,0 +1,146 @@ +# ============================================================================== +# pyvale: the python validation engine +# License: MIT +# Copyright (C) 2025 The Computer Aided Validation Team +# ============================================================================== + +""" +Specklegen: Speckle pattern generation using random disk placement with checking for overlap +================================================================================ +Script to generate a synthetic speckle pattern made from randomly placed circular +speckles (disks), run diagnostics on the generated image, and save both the image +and diagnostics to the selected folder. +""" + +import numpy as np +import time +import json +import os +import pyvale.specklegen as specklegen + +#%% +# Here we parse command line arguments to set the speckle pattern parameters. +# For ease of use in this example script we set parameter values directly in the +# code rather than via bash script. +# The parameter responsible for reducing overlap is set to 'True' in this example. +# Additionally, we have an extra parameter specifying the maximum number of attempts to place each speckle without overlap. + +speckle_size = 20 +screen_size_width = 1000 +screen_size_height = 800 +bit_depth = 8 +theme = 'white_on_black' +seed = 10 +sigma = 4.0 +reduce_overlap = True +type_gen = "random_disks" +attempts_tot = 300 +output_path = "src/pyvale/examples/specklegen/output/ex1b" + +print('Start') + +assert theme in ['black_on_white', 'white_on_black'], "Theme should be either 'black_on_white' or 'white_on_black'." + +if reduce_overlap: + print("Reducing overlap between speckles") +else: + print("Not reducing overlap between speckles") + +subfolder = f"/{type_gen}_{speckle_size}_{screen_size_width}_{screen_size_height}_{bit_depth}_{theme}_{seed}_{sigma}_{reduce_overlap}" +print(subfolder) +save_path = output_path + subfolder +if not os.path.exists(save_path): + os.makedirs(save_path) + +#%% +# We calculate parameteres aiming for approximately 50/50 black-to-white ratio. +# We now generate the speckle pattern using the specified parameters. +# The background and foreground colours are set based on the chosen theme and bit depth. +# We simply pass on one additional parameter to the function. + +speckle_area = np.pi * (speckle_size / 2) ** 2 +total_area = screen_size_width * screen_size_height +total_speckles = int((0.5 * total_area) / speckle_area) +print(f"Total number of speckles = {total_speckles}") +dynamic_range: int = 2**bit_depth - 1 +background_colour = 0 if theme == 'white_on_black' else dynamic_range +foreground_colour = dynamic_range if theme == 'white_on_black' else 0 + +feature_size_width = speckle_size +feature_size_height = speckle_size + +time_start = time.time() +image, results = specklegen.generate_speckles(screen_size_width, screen_size_height, + feature_size_width, feature_size_height, + foreground_colour, background_colour, + bit_depth, type_gen, seed, + total_speckles=total_speckles, + reduce_overlap=reduce_overlap, + sigma=sigma, attempts_tot=attempts_tot) +time_end = time.time() +time_taken = time_end - time_start +print(f"Time taken for speckle generation: {np.round(time_taken, 3)} seconds") + +np.savetxt(f"{save_path}/speckle_placement_results.csv", results, delimiter=",", + header="speckle_number, attempts, overlap(1/0/2), cent_x, cent_y", comments='', fmt=['%d', '%d', '%d', '%.3f', '%.3f']) + +#%% +# Now we run diagnostics on the generated speckle pattern and save the results. +# Finally, we print out the key statistics to the console. +# The plots are saved in the provided output folder. However, the diagnostic function outputs the matplotlib figures and axes, +# so the plot formatting could be changed from the default one used by the function. +# The speckles are placed in such a way as to reduce the overlap between them as much as possible. +# As a consequence of this, it can be seen that the black-to-white ratio in this example is considerably closer to unity than in the previous example (ex1a). +# This can also be supported by the irradiance value histogram. +# The number of 0 irradiance values, corresponding to black colour, and the 255 values, corresponding to white colour, became much more equal. + +print("") +print('Starting speckle pattern diagnostics...') +results = specklegen.speckle_pattern_statistics(image, dynamic_range) +plots = specklegen.speckle_pattern_plots(image, dynamic_range, save_path) + +with open(f"{save_path}/speckle_pattern_diagnostics.json", 'w') as f: + json.dump(results, f, indent=4) + +ratio = results.get("black_white_ratio", None) +mean_gradient = results.get("mean_intensity_gradient", None) +std_dev = results.get("std_dev_irradiance", None) +avg = results.get("avg_irradiance", None) +contrast = results.get("contrast", None) +entropy = results.get("shannon_entropy", None) +peak_to_mean = results.get("peak_to_mean_ratio", None) +skew = results.get("skewness", None) +kurt = results.get("kurtosis", None) +avg_speckle_size_fwhm = results.get("avg_speckle_size_fwhm", None) +avg_speckle_size_e2 = results.get("avg_speckle_size_e2", None) +H_fit_stats = results.get("H_fit_stats", None) +V_fit_stats = results.get("V_fit_stats", None) + +print("") +print("Speckle statistics:") + +print(f"Black/White ratio: {np.round(ratio, 3)}") +print(f"Mean intensity gradient: {np.round(mean_gradient, 3)}") +print(f"Standard deviation of irradiance values: {np.round(std_dev, 3)}") +print(f"Average irradiance value: {np.round(avg, 3)}") +print(f"Contrast (std/mean): {np.round(contrast, 3)}") +print(f"Skewness: {np.round(skew, 3)}") +print(f"Kurtosis: {np.round(kurt, 3)}") +print(f"Shannon entropy: {np.round(entropy, 3)}") +print(f"Peak to mean ratio: {np.round(peak_to_mean, 3)}") +print(f"Average speckle size (full width at half maximum): {np.round(avg_speckle_size_fwhm, 3)} pixels") +print(f"Average speckle size (1/e^2): {np.round(avg_speckle_size_e2, 3)} pixels") +print(f"R_squared: Horisontal fit: {np.round(H_fit_stats['R_squared'], 3)}, Vertical fit: {np.round(V_fit_stats['R_squared'], 3)}") + +#%% +# Finally, the relative errors beetween the specified speckle size and the speckle size approximated using cautocovariance are calculated. +error = np.abs(avg_speckle_size_fwhm - speckle_size) * 100 / speckle_size +print(f"Percentage error between requested speckle size and measured speckle size from FWHM: {np.round(error, 3)} %") +error = np.abs(avg_speckle_size_e2 - speckle_size) * 100 / speckle_size +print(f"Percentage error between requested speckle size and measured speckle size from 1/e^2: {np.round(error, 3)} %") +np.save(f"{save_path}/image.npy", image) +print("") +print('End :)') +print("") +print("") +print("") \ No newline at end of file diff --git a/src/pyvale/examples/specklegen/ex1c_random_disks_grid.py b/src/pyvale/examples/specklegen/ex1c_random_disks_grid.py new file mode 100644 index 00000000..c73da8ee --- /dev/null +++ b/src/pyvale/examples/specklegen/ex1c_random_disks_grid.py @@ -0,0 +1,146 @@ +# ============================================================================== +# pyvale: the python validation engine +# License: MIT +# Copyright (C) 2025 The Computer Aided Validation Team +# ============================================================================== + +""" +Specklegen: Speckle pattern generation using random disk-shaped speckle placement perturbation from a grid of regularly-placed disk-shaped speckles +================================================================================ +Script to generate a synthetic speckle pattern made from by randomly perturbating a grid of regularly-placed +disk-shaped speckles based on disrete uniform probability distribution, run diagnostics on the generated image, +and save both the image and diagnostics to the selected folder. +""" + +import numpy as np +import time +import json +import os +import pyvale.specklegen as specklegen + +#%% +# Here we parse command line arguments to set the speckle pattern parameters. +# For ease of use in this example script we set parameter values directly in the +# code rather than via bash script. +# The parameter responsible for reducing overlap is set to 'False' in this example. +# Here we select a different type of speckle generation compared with the previous examples (ex1a and ex1b). +# Additionally, we have an extra parameter specifying the maximum amount to move speckles by during grid perturbation (in pixels). + +speckle_size = 20 +screen_size_width = 1000 +screen_size_height = 800 +bit_depth = 8 +theme = 'white_on_black' +seed = 10 +sigma = 4.0 +reduce_overlap = False +perturbation_max = 12 +type_gen = "random_disks_grid" +output_path = "src/pyvale/examples/specklegen/output/ex1c" + +print('Start') + +assert theme in ['black_on_white', 'white_on_black'], "Theme should be either 'black_on_white' or 'white_on_black'." + +if reduce_overlap: + print("Reducing overlap between speckles") +else: + print("Not reducing overlap between speckles") + +subfolder = f"/{type_gen}_{speckle_size}_{screen_size_width}_{screen_size_height}_{bit_depth}_{theme}_{seed}" +print(subfolder) +save_path = output_path + subfolder +if not os.path.exists(save_path): + os.makedirs(save_path) + +#%% +# We calculate parameteres aiming for approximately 50/50 black-to-white ratio. +# We now generate the speckle pattern using the specified parameters. +# The background and foreground colours are set based on the chosen theme and bit depth. +# We simply pass on one additional parameter to the function. + +speckle_area = np.pi * (speckle_size / 2) ** 2 +total_area = screen_size_width * screen_size_height +total_speckles = int((0.5 * total_area) / speckle_area) +print(f"Total number of speckles = {total_speckles}") +dynamic_range: int = 2**bit_depth - 1 +background_colour = 0 if theme == 'white_on_black' else dynamic_range +foreground_colour = dynamic_range if theme == 'white_on_black' else 0 + +feature_size_width = speckle_size +feature_size_height = speckle_size + +time_start = time.time() +image, results = specklegen.generate_speckles(screen_size_width, screen_size_height, + feature_size_width, feature_size_height, + foreground_colour, background_colour, + bit_depth, type_gen, seed, + total_speckles=total_speckles, + reduce_overlap=reduce_overlap, + sigma=sigma, perturbation_max=perturbation_max) +time_end = time.time() +time_taken = time_end - time_start +print(f"Time taken for speckle generation: {np.round(time_taken, 3)} seconds") + +# save the speckle placement results +np.savetxt(f"{save_path}/speckle_placement_results.csv", results, delimiter=",", + header="speckle_number, attempts, overlap(1/0/2), cent_x, cent_y", comments='', fmt=['%d', '%d', '%d', '%.3f', '%.3f']) + +#%% +# Now we run diagnostics on the generated speckle pattern and save the results. +# Finally, we print out the key statistics to the console. +# The plots are saved in the provided output folder. However, the diagnostic function outputs the matplotlib figures and axes, +# so the plot formatting could be changed from the default one used by the function. +# The black-to-white ratio is better than achieved in ex1a, when we don’t check for speckle overlap. +# However, it is worse than the value achieved in ex1b, when we do check for speckle overlap. +# On the other hand, we still get the benefit of the improved black-to-white ratio at the reduced computational cost, +# as the runtime in this example is shorter than in ex1b. + +print("") +print('Starting speckle pattern diagnostics...') +results = specklegen.speckle_pattern_statistics(image, dynamic_range) +plots = specklegen.speckle_pattern_plots(image, dynamic_range, save_path) + +with open(f"{save_path}/speckle_pattern_diagnostics.json", 'w') as f: + json.dump(results, f, indent=4) + +ratio = results.get("black_white_ratio", None) +mean_gradient = results.get("mean_intensity_gradient", None) +std_dev = results.get("std_dev_irradiance", None) +avg = results.get("avg_irradiance", None) +contrast = results.get("contrast", None) +entropy = results.get("shannon_entropy", None) +peak_to_mean = results.get("peak_to_mean_ratio", None) +skew = results.get("skewness", None) +kurt = results.get("kurtosis", None) +avg_speckle_size_fwhm = results.get("avg_speckle_size_fwhm", None) +avg_speckle_size_e2 = results.get("avg_speckle_size_e2", None) +H_fit_stats = results.get("H_fit_stats", None) +V_fit_stats = results.get("V_fit_stats", None) + +print("") +print("Speckle statistics:") + +print(f"Black/White ratio: {np.round(ratio, 3)}") +print(f"Mean intensity gradient: {np.round(mean_gradient, 3)}") +print(f"Standard deviation of irradiance values: {np.round(std_dev, 3)}") +print(f"Average irradiance value: {np.round(avg, 3)}") +print(f"Contrast (std/mean): {np.round(contrast, 3)}") +print(f"Skewness: {np.round(skew, 3)}") +print(f"Kurtosis: {np.round(kurt, 3)}") +print(f"Shannon entropy: {np.round(entropy, 3)}") +print(f"Peak to mean ratio: {np.round(peak_to_mean, 3)}") +print(f"Average speckle size (full width at half maximum): {np.round(avg_speckle_size_fwhm, 3)} pixels") +print(f"Average speckle size (1/e^2): {np.round(avg_speckle_size_e2, 3)} pixels") +print(f"R_squared: Horisontal fit: {np.round(H_fit_stats['R_squared'], 3)}, Vertical fit: {np.round(V_fit_stats['R_squared'], 3)}") + +#%% +# Finally, the relative errors beetween the specified speckle size and the speckle size approximated using cautocovariance are calculated. +error = np.abs(avg_speckle_size_fwhm - speckle_size) * 100 / speckle_size +print(f"Percentage error between requested speckle size and measured speckle size from FWHM: {np.round(error, 3)} %") +error = np.abs(avg_speckle_size_e2 - speckle_size) * 100 / speckle_size +print(f"Percentage error between requested speckle size and measured speckle size from 1/e^2: {np.round(error, 3)} %") +np.save(f"{save_path}/image.npy", image) +print("") +print('End :)') +print("") \ No newline at end of file diff --git a/src/pyvale/examples/specklegen/ex2a_perlin_noise.py b/src/pyvale/examples/specklegen/ex2a_perlin_noise.py new file mode 100644 index 00000000..fbfb704d --- /dev/null +++ b/src/pyvale/examples/specklegen/ex2a_perlin_noise.py @@ -0,0 +1,129 @@ +# ============================================================================== +# pyvale: the python validation engine +# License: MIT +# Copyright (C) 2025 The Computer Aided Validation Team +# ============================================================================== + +""" +Specklegen: Speckle pattern generation using isotropic Perlin noise +================================================================================ +Script to generate a synthetic speckle pattern made using isotropic Perlin noise, run diagnostics on the generated image, and save both the image +and diagnostics to the selected folder. Isotroic Perlin noise in this case means that the speckle size is the same in both horisontal and vertical directions. + +This is a gradient-based noise algorithm that generates smooth and continuous random patterns. +It produces a texture with gradually occurring transitions. +Perlin noise achieves this by assigning random gradient vectors to grid points and then smoothly interpolating between them to create natural-looking transitions. +""" + +import numpy as np +import time +import json +import os +import pyvale.specklegen as specklegen + +#%% +# Here we parse command line arguments to set the speckle pattern parameters. +# For ease of use in this example script we set parameter values directly in the +# code rather than via bash script. +# The parameters are set to generate isotropic Perlin noise in this example. +# Perlin noise is defined by the number of noise periods width- and height-wise. +# They can be calculated from the corresponding width and height of a feature, which is a speckle size in our case, together with a screen size. +# The noise period is obtained by dividing a screen size by a speckle size. +# It should be noted that the screen size should be a multiple of the noise period number, otherwise the function would produce an error. + +speckle_size = 20 +screen_size_width = 1000 +screen_size_height = 800 +bit_depth = 8 +theme = 'white_on_black' +seed = 10 +type_gen = "perlin" +output_path = "src/pyvale/examples/specklegen/output/ex2a" + +print('Start') + +assert theme in ['black_on_white', 'white_on_black'], "Theme should be either 'black_on_white' or 'white_on_black'." + +subfolder = f"/{type_gen}_{speckle_size}_{screen_size_width}_{screen_size_height}_{bit_depth}_{theme}_{seed}" +print(subfolder) +save_path = output_path + subfolder +if not os.path.exists(save_path): + os.makedirs(save_path) + +#%% +# We now generate the speckle pattern using the specified parameters. +# The background and foreground colours are set based on the chosen theme and bit depth. +# It should be noted that there is no need to calculate the total number of speckles to generate like we did in the previous examples. + +dynamic_range: int = 2**bit_depth - 1 +background_colour = 0 if theme == 'white_on_black' else dynamic_range +foreground_colour = dynamic_range if theme == 'white_on_black' else 0 + +feature_size_width = speckle_size +feature_size_height = speckle_size + +time_start = time.time() +image = specklegen.generate_speckles(screen_size_width, screen_size_height, + feature_size_width, feature_size_height, + foreground_colour, background_colour, + bit_depth, type_gen, seed) +time_end = time.time() +time_taken = time_end - time_start +print(f"Time taken for speckle generation: {np.round(time_taken, 3)} seconds") + +#%% +# Now we run diagnostics on the generated speckle pattern and save the results. +# Finally, we print out the key statistics to the console. +# The plots are saved in the provided output folder. However, the diagnostic function outputs the matplotlib figures and axes, +# so the plot formatting could be changed from the default one used by the function. +# The black-to-white ratio is already close to unity, +# so there is no need to perform any additional operations related to speckle overlap reduction, like we did in ex1a, ex1b, and ex1c. + +print("") +print('Starting speckle pattern diagnostics...') +results = specklegen.speckle_pattern_statistics(image, dynamic_range) +plots = specklegen.speckle_pattern_plots(image, dynamic_range, save_path) + +with open(f"{save_path}/speckle_pattern_diagnostics.json", 'w') as f: + json.dump(results, f, indent=4) + +ratio = results.get("black_white_ratio", None) +mean_gradient = results.get("mean_intensity_gradient", None) +std_dev = results.get("std_dev_irradiance", None) +avg = results.get("avg_irradiance", None) +contrast = results.get("contrast", None) +entropy = results.get("shannon_entropy", None) +peak_to_mean = results.get("peak_to_mean_ratio", None) +skew = results.get("skewness", None) +kurt = results.get("kurtosis", None) +avg_speckle_size_fwhm = results.get("avg_speckle_size_fwhm", None) +avg_speckle_size_e2 = results.get("avg_speckle_size_e2", None) +H_fit_stats = results.get("H_fit_stats", None) +V_fit_stats = results.get("V_fit_stats", None) + +print("") +print("Speckle statistics:") + +print(f"Black/White ratio: {np.round(ratio, 3)}") +print(f"Mean intensity gradient: {np.round(mean_gradient, 3)}") +print(f"Standard deviation of irradiance values: {np.round(std_dev, 3)}") +print(f"Average irradiance value: {np.round(avg, 3)}") +print(f"Contrast (std/mean): {np.round(contrast, 3)}") +print(f"Skewness: {np.round(skew, 3)}") +print(f"Kurtosis: {np.round(kurt, 3)}") +print(f"Shannon entropy: {np.round(entropy, 3)}") +print(f"Peak to mean ratio: {np.round(peak_to_mean, 3)}") +print(f"Average speckle size (full width at half maximum): {np.round(avg_speckle_size_fwhm, 3)} pixels") +print(f"Average speckle size (1/e^2): {np.round(avg_speckle_size_e2, 3)} pixels") +print(f"R_squared: Horisontal fit: {np.round(H_fit_stats['R_squared'], 3)}, Vertical fit: {np.round(V_fit_stats['R_squared'], 3)}") + +#%% +# Finally, the relative errors beetween the specified speckle size and the speckle size approximated using cautocovariance are calculated. +error = np.abs(avg_speckle_size_fwhm - speckle_size) * 100 / speckle_size +print(f"Percentage error between requested speckle size and measured speckle size from FWHM: {np.round(error, 3)} %") +error = np.abs(avg_speckle_size_e2 - speckle_size) * 100 / speckle_size +print(f"Percentage error between requested speckle size and measured speckle size from 1/e^2: {np.round(error, 3)} %") +np.save(f"{save_path}/image.npy", image) +print("") +print('End :)') +print("") \ No newline at end of file diff --git a/src/pyvale/examples/specklegen/ex3a_fractal_noise.py b/src/pyvale/examples/specklegen/ex3a_fractal_noise.py new file mode 100644 index 00000000..e3e551a9 --- /dev/null +++ b/src/pyvale/examples/specklegen/ex3a_fractal_noise.py @@ -0,0 +1,139 @@ +# ============================================================================== +# pyvale: the python validation engine +# License: MIT +# Copyright (C) 2025 The Computer Aided Validation Team +# ============================================================================== + +""" +Specklegen: Speckle pattern generation using isotropic fractal noise +================================================================================ +Script to generate a synthetic speckle pattern made using isotropic fractal noise, run diagnostics on the generated image, and save both the image +and diagnostics to the selected folder. Isotroic fractal noise in this case means that the speckle size is the same in both horisontal and vertical directions. + +The fractal noise pattern is usually characterised by self-similarity across multiple scales. +It is commonly used to produce realistic, organic textures. +Fractal noise combines smooth but irregular variations at different levels of detail by layering several frequencies of Perlin noise (octaves), +each with its own amplitude. +""" + +import numpy as np +import time +import json +import os +import pyvale.specklegen as specklegen + +#%% +# Here we parse command line arguments to set the speckle pattern parameters. +# For ease of use in this example script we set parameter values directly in the +# code rather than via bash script. +# The parameters are set to generate isotropic fractal noise in this example. +# Fractal noise is defined by the number of noise periods width- and height-wise. +# They can be calculated from the corresponding width and height of a feature, which is a speckle size in our case, together with a screen size. +# The noise period (res) is obtained by dividing a screen size by a speckle size. +# Additionally, we see two new parameters here: octaves and lacunarity. +# Octaves define a number of detail levels in the generated pattern. +# Lacunarity is a frequency factor between two octaves. It essentially defines the amount of detail added or removed at each octave. +# It should be noted that the screen size should be a multiple of the lacunarity^(octaves-1)*res, otherwise the function would produce an error. +# Consequently, the restrictions placed on the relationship between octaves, lacunarity, speckle and screen sizes make the parameter definition more complex. + +speckle_size = 20 +screen_size_width = 1000 +screen_size_height = 800 +bit_depth = 8 +theme = 'white_on_black' +seed = 10 +type_gen = "fractal" +octaves = 3 +lacunarity = 2 +output_path = "src/pyvale/examples/specklegen/output/ex3a" + +print('Start') + +assert theme in ['black_on_white', 'white_on_black'], "Theme should be either 'black_on_white' or 'white_on_black'." + +subfolder = f"/{type_gen}_{speckle_size}_{screen_size_width}_{screen_size_height}_{bit_depth}_{theme}_{seed}" +print(subfolder) +save_path = output_path + subfolder +if not os.path.exists(save_path): + os.makedirs(save_path) + +#%% +# We now generate the speckle pattern using the specified parameters. +# The background and foreground colours are set based on the chosen theme and bit depth. +# We simply pass on the two additional parameters to the function. + +dynamic_range: int = 2**bit_depth - 1 +background_colour = 0 if theme == 'white_on_black' else dynamic_range +foreground_colour = dynamic_range if theme == 'white_on_black' else 0 + +feature_size_width = speckle_size +feature_size_height = speckle_size + +time_start = time.time() +image = specklegen.generate_speckles(screen_size_width, screen_size_height, + feature_size_width, feature_size_height, + foreground_colour, background_colour, + bit_depth, type_gen, seed, + octaves=octaves, lacunarity=lacunarity) +time_end = time.time() +time_taken = time_end - time_start +print(f"Time taken for speckle generation: {np.round(time_taken, 3)} seconds") + +#%% +# Now we run diagnostics on the generated speckle pattern and save the results. +# Finally, we print out the key statistics to the console. +# The plots are saved in the provided output folder. However, the diagnostic function outputs the matplotlib figures and axes, +# so the plot formatting could be changed from the default one used by the function. +# The black-to-white ratio is already close to unity, +# so there is no need to perform any additional operations related to speckle overlap reduction, like we did in ex1a, ex1b, and ex1c. +# Compared with the previous example (ex2a), the generated speckle pattern does exhibit more textured appearance. +# However, this makes the pattern less realistic and hence less suitable for our applications. + +print("") +print('Starting speckle pattern diagnostics...') +results = specklegen.speckle_pattern_statistics(image, dynamic_range) +plots = specklegen.speckle_pattern_plots(image, dynamic_range, save_path) + +with open(f"{save_path}/speckle_pattern_diagnostics.json", 'w') as f: + json.dump(results, f, indent=4) + +ratio = results.get("black_white_ratio", None) +mean_gradient = results.get("mean_intensity_gradient", None) +std_dev = results.get("std_dev_irradiance", None) +avg = results.get("avg_irradiance", None) +contrast = results.get("contrast", None) +entropy = results.get("shannon_entropy", None) +peak_to_mean = results.get("peak_to_mean_ratio", None) +skew = results.get("skewness", None) +kurt = results.get("kurtosis", None) +avg_speckle_size_fwhm = results.get("avg_speckle_size_fwhm", None) +avg_speckle_size_e2 = results.get("avg_speckle_size_e2", None) +H_fit_stats = results.get("H_fit_stats", None) +V_fit_stats = results.get("V_fit_stats", None) + +print("") +print("Speckle statistics:") + +print(f"Black/White ratio: {np.round(ratio, 3)}") +print(f"Mean intensity gradient: {np.round(mean_gradient, 3)}") +print(f"Standard deviation of irradiance values: {np.round(std_dev, 3)}") +print(f"Average irradiance value: {np.round(avg, 3)}") +print(f"Contrast (std/mean): {np.round(contrast, 3)}") +print(f"Skewness: {np.round(skew, 3)}") +print(f"Kurtosis: {np.round(kurt, 3)}") +print(f"Shannon entropy: {np.round(entropy, 3)}") +print(f"Peak to mean ratio: {np.round(peak_to_mean, 3)}") +print(f"Average speckle size (full width at half maximum): {np.round(avg_speckle_size_fwhm, 3)} pixels") +print(f"Average speckle size (1/e^2): {np.round(avg_speckle_size_e2, 3)} pixels") +print(f"R_squared: Horisontal fit: {np.round(H_fit_stats['R_squared'], 3)}, Vertical fit: {np.round(V_fit_stats['R_squared'], 3)}") + +#%% +# Finally, the relative errors beetween the specified speckle size and the speckle size approximated using cautocovariance are calculated. +error = np.abs(avg_speckle_size_fwhm - speckle_size) * 100 / speckle_size +print(f"Percentage error between requested speckle size and measured speckle size from FWHM: {np.round(error, 3)} %") +error = np.abs(avg_speckle_size_e2 - speckle_size) * 100 / speckle_size +print(f"Percentage error between requested speckle size and measured speckle size from 1/e^2: {np.round(error, 3)} %") +np.save(f"{save_path}/image.npy", image) +print("") +print('End :)') +print("") \ No newline at end of file diff --git a/src/pyvale/examples/specklegen/ex4a_simplex_noise.py b/src/pyvale/examples/specklegen/ex4a_simplex_noise.py new file mode 100644 index 00000000..16d430da --- /dev/null +++ b/src/pyvale/examples/specklegen/ex4a_simplex_noise.py @@ -0,0 +1,123 @@ +# ============================================================================== +# pyvale: the python validation engine +# License: MIT +# Copyright (C) 2025 The Computer Aided Validation Team +# ============================================================================== + +""" +Specklegen: Speckle pattern generation using isotropic Simplex noise +================================================================================ +Script to generate a synthetic speckle pattern made using isotropic Simplex noise, run diagnostics on the generated image, and save both the image +and diagnostics to the selected folder. Isotroic Simplex noise in this case means that the speckle size is the same in both horisontal and vertical directions. + +Simplex noise is an enhanced version of Perlin noise that aims to produce more consistent and isotropic noise patterns. +""" + +import numpy as np +import time +import json +import os +import pyvale.specklegen as specklegen + +#%% +# Here we parse command line arguments to set the speckle pattern parameters. +# For ease of use in this example script we set parameter values directly in the +# code rather than via bash script. +# The parameters are set to generate isotropic Simplex noise in this example. +# Contrary to the Perlin (ex2a) and fractal (ex3a) noise, there are no restrictions placed on the parameters defining Simplex noise. + +speckle_size = 20 +screen_size_width = 1000 +screen_size_height = 800 +bit_depth = 8 +theme = 'white_on_black' +seed = 10 +type_gen = "simplex" +output_path = "src/pyvale/examples/specklegen/output/ex4a" + +print('Start') + +assert theme in ['black_on_white', 'white_on_black'], "Theme should be either 'black_on_white' or 'white_on_black'." + +subfolder = f"/{type_gen}_{speckle_size}_{screen_size_width}_{screen_size_height}_{bit_depth}_{theme}_{seed}" +print(subfolder) +save_path = output_path + subfolder +if not os.path.exists(save_path): + os.makedirs(save_path) + +#%% +# We now generate the speckle pattern using the specified parameters. +# The background and foreground colours are set based on the chosen theme and bit depth. + +dynamic_range: int = 2**bit_depth - 1 +background_colour = 0 if theme == 'white_on_black' else dynamic_range +foreground_colour = dynamic_range if theme == 'white_on_black' else 0 + +feature_size_width = speckle_size +feature_size_height = speckle_size + +time_start = time.time() +image = specklegen.generate_speckles(screen_size_width, screen_size_height, + feature_size_width, feature_size_height, + foreground_colour, background_colour, + bit_depth, type_gen, seed) +time_end = time.time() +time_taken = time_end - time_start +print(f"Time taken for speckle generation: {np.round(time_taken, 3)} seconds") + +#%% +# Now we run diagnostics on the generated speckle pattern and save the results. +# Finally, we print out the key statistics to the console. +# The plots are saved in the provided output folder. However, the diagnostic function outputs the matplotlib figures and axes, +# so the plot formatting could be changed from the default one used by the function. +# The black-to-white ratio is already close to unity, +# so there is no need to perform any additional operations related to speckle overlap reduction, like we did in ex1a, ex1b, and ex1c. + +print("") +print('Starting speckle pattern diagnostics...') +results = specklegen.speckle_pattern_statistics(image, dynamic_range) +plots = specklegen.speckle_pattern_plots(image, dynamic_range, save_path) + +with open(f"{save_path}/speckle_pattern_diagnostics.json", 'w') as f: + json.dump(results, f, indent=4) + +ratio = results.get("black_white_ratio", None) +mean_gradient = results.get("mean_intensity_gradient", None) +std_dev = results.get("std_dev_irradiance", None) +avg = results.get("avg_irradiance", None) +contrast = results.get("contrast", None) +entropy = results.get("shannon_entropy", None) +peak_to_mean = results.get("peak_to_mean_ratio", None) +skew = results.get("skewness", None) +kurt = results.get("kurtosis", None) +avg_speckle_size_fwhm = results.get("avg_speckle_size_fwhm", None) +avg_speckle_size_e2 = results.get("avg_speckle_size_e2", None) +H_fit_stats = results.get("H_fit_stats", None) +V_fit_stats = results.get("V_fit_stats", None) + +print("") +print("Speckle statistics:") + +print(f"Black/White ratio: {np.round(ratio, 3)}") +print(f"Mean intensity gradient: {np.round(mean_gradient, 3)}") +print(f"Standard deviation of irradiance values: {np.round(std_dev, 3)}") +print(f"Average irradiance value: {np.round(avg, 3)}") +print(f"Contrast (std/mean): {np.round(contrast, 3)}") +print(f"Skewness: {np.round(skew, 3)}") +print(f"Kurtosis: {np.round(kurt, 3)}") +print(f"Shannon entropy: {np.round(entropy, 3)}") +print(f"Peak to mean ratio: {np.round(peak_to_mean, 3)}") +print(f"Average speckle size (full width at half maximum): {np.round(avg_speckle_size_fwhm, 3)} pixels") +print(f"Average speckle size (1/e^2): {np.round(avg_speckle_size_e2, 3)} pixels") +print(f"R_squared: Horisontal fit: {np.round(H_fit_stats['R_squared'], 3)}, Vertical fit: {np.round(V_fit_stats['R_squared'], 3)}") + +#%% +# Finally, the relative errors beetween the specified speckle size and the speckle size approximated using cautocovariance are calculated. +error = np.abs(avg_speckle_size_fwhm - speckle_size) * 100 / speckle_size +print(f"Percentage error between requested speckle size and measured speckle size from FWHM: {np.round(error, 3)} %") +error = np.abs(avg_speckle_size_e2 - speckle_size) * 100 / speckle_size +print(f"Percentage error between requested speckle size and measured speckle size from 1/e^2: {np.round(error, 3)} %") +np.save(f"{save_path}/image.npy", image) +print("") +print('End :)') +print("") \ No newline at end of file diff --git a/src/pyvale/examples/specklegen/output/ex1a/random_disks_20_1000_800_8_white_on_black_10_4.0_False/autocovariance.jpg b/src/pyvale/examples/specklegen/output/ex1a/random_disks_20_1000_800_8_white_on_black_10_4.0_False/autocovariance.jpg new file mode 100644 index 00000000..58fda662 Binary files /dev/null and b/src/pyvale/examples/specklegen/output/ex1a/random_disks_20_1000_800_8_white_on_black_10_4.0_False/autocovariance.jpg differ diff --git a/src/pyvale/examples/specklegen/output/ex1a/random_disks_20_1000_800_8_white_on_black_10_4.0_False/frequency_spectrum.jpg b/src/pyvale/examples/specklegen/output/ex1a/random_disks_20_1000_800_8_white_on_black_10_4.0_False/frequency_spectrum.jpg new file mode 100644 index 00000000..33a31c3f Binary files /dev/null and b/src/pyvale/examples/specklegen/output/ex1a/random_disks_20_1000_800_8_white_on_black_10_4.0_False/frequency_spectrum.jpg differ diff --git a/src/pyvale/examples/specklegen/output/ex1a/random_disks_20_1000_800_8_white_on_black_10_4.0_False/pixel_value_histogram.jpg b/src/pyvale/examples/specklegen/output/ex1a/random_disks_20_1000_800_8_white_on_black_10_4.0_False/pixel_value_histogram.jpg new file mode 100644 index 00000000..b98cb74c Binary files /dev/null and b/src/pyvale/examples/specklegen/output/ex1a/random_disks_20_1000_800_8_white_on_black_10_4.0_False/pixel_value_histogram.jpg differ diff --git a/src/pyvale/examples/specklegen/output/ex1a/random_disks_20_1000_800_8_white_on_black_10_4.0_False/speckle_pattern.jpg b/src/pyvale/examples/specklegen/output/ex1a/random_disks_20_1000_800_8_white_on_black_10_4.0_False/speckle_pattern.jpg new file mode 100644 index 00000000..58d9c12a Binary files /dev/null and b/src/pyvale/examples/specklegen/output/ex1a/random_disks_20_1000_800_8_white_on_black_10_4.0_False/speckle_pattern.jpg differ diff --git a/src/pyvale/examples/specklegen/output/ex1b/random_disks_20_1000_800_8_white_on_black_10_4.0_True/autocovariance.jpg b/src/pyvale/examples/specklegen/output/ex1b/random_disks_20_1000_800_8_white_on_black_10_4.0_True/autocovariance.jpg new file mode 100644 index 00000000..bbc17b56 Binary files /dev/null and b/src/pyvale/examples/specklegen/output/ex1b/random_disks_20_1000_800_8_white_on_black_10_4.0_True/autocovariance.jpg differ diff --git a/src/pyvale/examples/specklegen/output/ex1b/random_disks_20_1000_800_8_white_on_black_10_4.0_True/frequency_spectrum.jpg b/src/pyvale/examples/specklegen/output/ex1b/random_disks_20_1000_800_8_white_on_black_10_4.0_True/frequency_spectrum.jpg new file mode 100644 index 00000000..7871103e Binary files /dev/null and b/src/pyvale/examples/specklegen/output/ex1b/random_disks_20_1000_800_8_white_on_black_10_4.0_True/frequency_spectrum.jpg differ diff --git a/src/pyvale/examples/specklegen/output/ex1b/random_disks_20_1000_800_8_white_on_black_10_4.0_True/pixel_value_histogram.jpg b/src/pyvale/examples/specklegen/output/ex1b/random_disks_20_1000_800_8_white_on_black_10_4.0_True/pixel_value_histogram.jpg new file mode 100644 index 00000000..4fa73691 Binary files /dev/null and b/src/pyvale/examples/specklegen/output/ex1b/random_disks_20_1000_800_8_white_on_black_10_4.0_True/pixel_value_histogram.jpg differ diff --git a/src/pyvale/examples/specklegen/output/ex1b/random_disks_20_1000_800_8_white_on_black_10_4.0_True/speckle_pattern.jpg b/src/pyvale/examples/specklegen/output/ex1b/random_disks_20_1000_800_8_white_on_black_10_4.0_True/speckle_pattern.jpg new file mode 100644 index 00000000..c8c15813 Binary files /dev/null and b/src/pyvale/examples/specklegen/output/ex1b/random_disks_20_1000_800_8_white_on_black_10_4.0_True/speckle_pattern.jpg differ diff --git a/src/pyvale/examples/specklegen/output/ex1c/random_disks_grid_20_1000_800_8_white_on_black_10/autocovariance.jpg b/src/pyvale/examples/specklegen/output/ex1c/random_disks_grid_20_1000_800_8_white_on_black_10/autocovariance.jpg new file mode 100644 index 00000000..0bde577b Binary files /dev/null and b/src/pyvale/examples/specklegen/output/ex1c/random_disks_grid_20_1000_800_8_white_on_black_10/autocovariance.jpg differ diff --git a/src/pyvale/examples/specklegen/output/ex1c/random_disks_grid_20_1000_800_8_white_on_black_10/frequency_spectrum.jpg b/src/pyvale/examples/specklegen/output/ex1c/random_disks_grid_20_1000_800_8_white_on_black_10/frequency_spectrum.jpg new file mode 100644 index 00000000..982f8164 Binary files /dev/null and b/src/pyvale/examples/specklegen/output/ex1c/random_disks_grid_20_1000_800_8_white_on_black_10/frequency_spectrum.jpg differ diff --git a/src/pyvale/examples/specklegen/output/ex1c/random_disks_grid_20_1000_800_8_white_on_black_10/pixel_value_histogram.jpg b/src/pyvale/examples/specklegen/output/ex1c/random_disks_grid_20_1000_800_8_white_on_black_10/pixel_value_histogram.jpg new file mode 100644 index 00000000..f858752c Binary files /dev/null and b/src/pyvale/examples/specklegen/output/ex1c/random_disks_grid_20_1000_800_8_white_on_black_10/pixel_value_histogram.jpg differ diff --git a/src/pyvale/examples/specklegen/output/ex1c/random_disks_grid_20_1000_800_8_white_on_black_10/speckle_pattern.jpg b/src/pyvale/examples/specklegen/output/ex1c/random_disks_grid_20_1000_800_8_white_on_black_10/speckle_pattern.jpg new file mode 100644 index 00000000..9ce6989e Binary files /dev/null and b/src/pyvale/examples/specklegen/output/ex1c/random_disks_grid_20_1000_800_8_white_on_black_10/speckle_pattern.jpg differ diff --git a/src/pyvale/examples/specklegen/output/ex2a/perlin_20_1000_800_8_white_on_black_10/autocovariance.jpg b/src/pyvale/examples/specklegen/output/ex2a/perlin_20_1000_800_8_white_on_black_10/autocovariance.jpg new file mode 100644 index 00000000..272be0ce Binary files /dev/null and b/src/pyvale/examples/specklegen/output/ex2a/perlin_20_1000_800_8_white_on_black_10/autocovariance.jpg differ diff --git a/src/pyvale/examples/specklegen/output/ex2a/perlin_20_1000_800_8_white_on_black_10/frequency_spectrum.jpg b/src/pyvale/examples/specklegen/output/ex2a/perlin_20_1000_800_8_white_on_black_10/frequency_spectrum.jpg new file mode 100644 index 00000000..0882cacd Binary files /dev/null and b/src/pyvale/examples/specklegen/output/ex2a/perlin_20_1000_800_8_white_on_black_10/frequency_spectrum.jpg differ diff --git a/src/pyvale/examples/specklegen/output/ex2a/perlin_20_1000_800_8_white_on_black_10/pixel_value_histogram.jpg b/src/pyvale/examples/specklegen/output/ex2a/perlin_20_1000_800_8_white_on_black_10/pixel_value_histogram.jpg new file mode 100644 index 00000000..8bd9ea35 Binary files /dev/null and b/src/pyvale/examples/specklegen/output/ex2a/perlin_20_1000_800_8_white_on_black_10/pixel_value_histogram.jpg differ diff --git a/src/pyvale/examples/specklegen/output/ex2a/perlin_20_1000_800_8_white_on_black_10/speckle_pattern.jpg b/src/pyvale/examples/specklegen/output/ex2a/perlin_20_1000_800_8_white_on_black_10/speckle_pattern.jpg new file mode 100644 index 00000000..acc267cc Binary files /dev/null and b/src/pyvale/examples/specklegen/output/ex2a/perlin_20_1000_800_8_white_on_black_10/speckle_pattern.jpg differ diff --git a/src/pyvale/examples/specklegen/output/ex3a/fractal_20_1000_800_8_white_on_black_10/autocovariance.jpg b/src/pyvale/examples/specklegen/output/ex3a/fractal_20_1000_800_8_white_on_black_10/autocovariance.jpg new file mode 100644 index 00000000..c06754d9 Binary files /dev/null and b/src/pyvale/examples/specklegen/output/ex3a/fractal_20_1000_800_8_white_on_black_10/autocovariance.jpg differ diff --git a/src/pyvale/examples/specklegen/output/ex3a/fractal_20_1000_800_8_white_on_black_10/frequency_spectrum.jpg b/src/pyvale/examples/specklegen/output/ex3a/fractal_20_1000_800_8_white_on_black_10/frequency_spectrum.jpg new file mode 100644 index 00000000..86ec7d64 Binary files /dev/null and b/src/pyvale/examples/specklegen/output/ex3a/fractal_20_1000_800_8_white_on_black_10/frequency_spectrum.jpg differ diff --git a/src/pyvale/examples/specklegen/output/ex3a/fractal_20_1000_800_8_white_on_black_10/pixel_value_histogram.jpg b/src/pyvale/examples/specklegen/output/ex3a/fractal_20_1000_800_8_white_on_black_10/pixel_value_histogram.jpg new file mode 100644 index 00000000..29e52d12 Binary files /dev/null and b/src/pyvale/examples/specklegen/output/ex3a/fractal_20_1000_800_8_white_on_black_10/pixel_value_histogram.jpg differ diff --git a/src/pyvale/examples/specklegen/output/ex3a/fractal_20_1000_800_8_white_on_black_10/speckle_pattern.jpg b/src/pyvale/examples/specklegen/output/ex3a/fractal_20_1000_800_8_white_on_black_10/speckle_pattern.jpg new file mode 100644 index 00000000..0e349505 Binary files /dev/null and b/src/pyvale/examples/specklegen/output/ex3a/fractal_20_1000_800_8_white_on_black_10/speckle_pattern.jpg differ diff --git a/src/pyvale/examples/specklegen/output/ex4a/simplex_20_1000_800_8_white_on_black_10/autocovariance.jpg b/src/pyvale/examples/specklegen/output/ex4a/simplex_20_1000_800_8_white_on_black_10/autocovariance.jpg new file mode 100644 index 00000000..04c1f2ae Binary files /dev/null and b/src/pyvale/examples/specklegen/output/ex4a/simplex_20_1000_800_8_white_on_black_10/autocovariance.jpg differ diff --git a/src/pyvale/examples/specklegen/output/ex4a/simplex_20_1000_800_8_white_on_black_10/frequency_spectrum.jpg b/src/pyvale/examples/specklegen/output/ex4a/simplex_20_1000_800_8_white_on_black_10/frequency_spectrum.jpg new file mode 100644 index 00000000..c6c840cb Binary files /dev/null and b/src/pyvale/examples/specklegen/output/ex4a/simplex_20_1000_800_8_white_on_black_10/frequency_spectrum.jpg differ diff --git a/src/pyvale/examples/specklegen/output/ex4a/simplex_20_1000_800_8_white_on_black_10/pixel_value_histogram.jpg b/src/pyvale/examples/specklegen/output/ex4a/simplex_20_1000_800_8_white_on_black_10/pixel_value_histogram.jpg new file mode 100644 index 00000000..62bbf7f3 Binary files /dev/null and b/src/pyvale/examples/specklegen/output/ex4a/simplex_20_1000_800_8_white_on_black_10/pixel_value_histogram.jpg differ diff --git a/src/pyvale/examples/specklegen/output/ex4a/simplex_20_1000_800_8_white_on_black_10/speckle_pattern.jpg b/src/pyvale/examples/specklegen/output/ex4a/simplex_20_1000_800_8_white_on_black_10/speckle_pattern.jpg new file mode 100644 index 00000000..128f89de Binary files /dev/null and b/src/pyvale/examples/specklegen/output/ex4a/simplex_20_1000_800_8_white_on_black_10/speckle_pattern.jpg differ diff --git a/src/pyvale/sensorsim/cython/rastercyth.c b/src/pyvale/sensorsim/cython/rastercyth.c index 6679b499..c60e916c 100644 --- a/src/pyvale/sensorsim/cython/rastercyth.c +++ b/src/pyvale/sensorsim/cython/rastercyth.c @@ -418,6 +418,9 @@ END: Cython Metadata */ enum { __pyx_check_sizeof_voidp = 1 / (int)(SIZEOF_VOID_P == sizeof(void*)) }; #endif #endif +#ifndef CYTHON_LOCK_AND_GIL_DEADLOCK_AVOIDANCE_TIME + #define CYTHON_LOCK_AND_GIL_DEADLOCK_AVOIDANCE_TIME 100 +#endif #ifndef __has_attribute #define __has_attribute(x) 0 #endif @@ -1298,6 +1301,7 @@ static CYTHON_INLINE Py_hash_t __Pyx_PyIndex_AsHash_t(PyObject*); typedef sdigit __Pyx_compact_pylong; typedef digit __Pyx_compact_upylong; #endif + static CYTHON_INLINE int __Pyx_PyLong_CompactAsLong(PyObject *x, long *return_value); #if PY_VERSION_HEX >= 0x030C00A5 #define __Pyx_PyLong_Digits(x) (((PyLongObject*)x)->long_value.ob_digit) #else @@ -32076,6 +32080,17 @@ static CYTHON_INLINE PyObject * __Pyx_PyBool_FromLong(long b) { static CYTHON_INLINE PyObject * __Pyx_PyLong_FromSize_t(size_t ival) { return PyLong_FromSize_t(ival); } +#if CYTHON_USE_PYLONG_INTERNALS +static CYTHON_INLINE int __Pyx_PyLong_CompactAsLong(PyObject *x, long *return_value) { + if (unlikely(!__Pyx_PyLong_IsCompact(x))) + return 0; + Py_ssize_t value = __Pyx_PyLong_CompactValue(x); + if ((sizeof(long) < sizeof(Py_ssize_t)) && unlikely(value != (long) value)) + return 0; + *return_value = (long) value; + return 1; +} +#endif /* MultiPhaseInitModuleState */ diff --git a/src/pyvale/specklegen/__init__.py b/src/pyvale/specklegen/__init__.py new file mode 100644 index 00000000..d55ebac2 --- /dev/null +++ b/src/pyvale/specklegen/__init__.py @@ -0,0 +1,8 @@ +#=============================================================================== +# pyvale: the python validation engine +# License: MIT +# Copyright (C) 2025 The Computer Aided Validation Team +#=============================================================================== + +from .speckle_generation import * +from .speckle_pattern_diagnostics import * \ No newline at end of file diff --git a/src/pyvale/specklegen/speckle_generation.py b/src/pyvale/specklegen/speckle_generation.py new file mode 100644 index 00000000..5c81f9f2 --- /dev/null +++ b/src/pyvale/specklegen/speckle_generation.py @@ -0,0 +1,489 @@ +import numpy as np +from scipy import ndimage +from perlin_numpy import ( + generate_perlin_noise_2d, generate_fractal_noise_2d) +import opensimplex as simplex + +def pixelsInDisk(cent_x: int, cent_y: int, + screen_size_width: int, screen_size_height: int, + speckle_size: float, foreground_colour: int, image: np.ndarray, + check_overlap: bool) -> int: + """A function to set pixels in a disk shape on the image array + + Parameters + ---------- + cent_x : int + Center x coordinate of the disk + cent_y : int + Center y coordinate of the disk + screen_size_width : int + Dimension of the image width-wise (pixels) + screen_size_height : int + Dimension of the image height-wise (pixels) + speckle_size : float + Diameter of the speckle disk (pixels) + foreground_colour : int + Colour value to set for the speckle pixels + image : np.ndarray + 2D numpy array representing the image + check_overlap : bool + If True, check for overlap with existing foreground pixels and does not modify the image + If False, sets the pixels in the disk shape to the foreground colour + + Returns + ------- + int + 0 if overlap detected, 1 if no overlap (when check_overlap is True), otherwise a function modifies the image + """ + + box_max_x = int(np.ceil(cent_x + speckle_size / 2)) + box_min_x = int(np.floor(cent_x - speckle_size / 2)) + + box_max_y = int(np.ceil(cent_y + speckle_size / 2)) + box_min_y = int(np.floor(cent_y - speckle_size / 2)) + + if check_overlap: + over_lap = False + + for yy in range(box_min_y, box_max_y): + for xx in range(box_min_x, box_max_x): + if yy < 0 or yy >= screen_size_height or xx < 0 or xx >= screen_size_width: + continue + distance = (xx + 0.5 - cent_x)**2 + (yy + 0.5 - cent_y)**2 + if distance <= (speckle_size / 2)**2: + if image[yy, xx] == foreground_colour: + over_lap = True + break + if over_lap: + # print(f"Overlap detected at ({xx}, {yy})") + break + return 0 if over_lap else 1 + + else: + for yy in range(box_min_y, box_max_y): + for xx in range(box_min_x, box_max_x): + if yy < 0 or yy >= screen_size_height or xx < 0 or xx >= screen_size_width: + continue + distance = (xx + 0.5 - cent_x)**2 + (yy + 0.5 - cent_y)**2 + if distance <= (speckle_size / 2)**2: + image[yy, xx] = foreground_colour + + +def generate_speckles_random_disks(screen_size_width: int, screen_size_height: int, + feature_size_width: float, feature_size_height: float, + foreground_colour: int, background_colour: int, + bit_depth: int, type_gen: str, seed: int, **kwargs) -> np.ndarray: + """A function to generate a speckle pattern image by randomly placing disk-shaped speckles + based on uniform probability distribution. + + Parameters + ---------- + screen_size_width : int + Dimension of the image width-wise (pixels) + screen_size_height : int + Dimension of the image height-wise (pixels) + feature_size_width : float + Speckle size width-wise + feature_size_height : float + Speckle size height-wise + foreground_colour : int + Colour value for the speckles + background_colour : int + Colour value for the background + bit_depth : int + Bit depth of the image (8 or 16) + type_gen : str + Type of noise to generate (must be 'random_disks' for this function) + seed : int + Random seed for the noise generation + total_speckles : int + Total number of speckles to place + reduce_overlap : bool + If True, attempts to reduce overlap between speckles + sigma : float + Standard deviation for Gaussian blur applied after speckle placement + attempts_tot : int (optional) + Maximum number of attempts to place each speckle without overlap (if reduce_overlap is True) + + Returns + ------- + np.ndarray + Speckle pattern image as a 2D numpy array and speckle generation stats + """ + + speckle_size = (feature_size_width + feature_size_height) / 2 + total_speckles = kwargs.get("total_speckles") + reduce_overlap = kwargs.get("reduce_overlap") + sigma = kwargs.get("sigma") + attempts_tot = kwargs.get("attempts_tot", 100) + + # check that non of the required parameters are None + assert speckle_size is not None, "speckle_size parameter is required." + assert total_speckles is not None, "total_speckles parameter is required." + assert reduce_overlap is not None, "reduce_overlap parameter is required." + assert sigma is not None, "sigma parameter is required." + assert type_gen == "random_disks", "Type_gen must be 'random_disks' for this function." + + np.random.seed(seed) + image: np.ndarray = np.ones((screen_size_height,screen_size_width), + dtype=np.uint16 if bit_depth == 16 else np.uint8) * background_colour + + results = np.zeros([total_speckles, 5]) # speckle number, attempts, with/without/not checked overlap (1/0/2), cent_x, cent_y + + if reduce_overlap: + for i in range(int(total_speckles)): + # reduce overlap by checking if the new speckle overlaps with existing ones + check_overlap = True + over_lap = True + attempts = 0 + # print(f"Placing speckle {i+1}/{total_speckles}") + while over_lap and attempts < attempts_tot: + # print(f"Attempt {attempts+1} to place speckle {i+1}") + over_lap = False + cent_x = np.random.uniform(0, screen_size_width) + cent_y = np.random.uniform(0, screen_size_height) + + # check colour + if image[int(cent_y), int(cent_x)] == foreground_colour: + over_lap = True + attempts += 1 + # print(f"Overlap detected at centre ({cent_x}, {cent_y})") + continue + over_lap_int = pixelsInDisk(cent_x, cent_y, + screen_size_width, screen_size_height, + speckle_size, foreground_colour, image, check_overlap) + over_lap = True if over_lap_int == 0 else False + + results[i, 0] = i + 1 + results[i, 1] = attempts + 1 + results[i, 2] = 0 if over_lap else 1 + results[i, 3] = cent_x + results[i, 4] = cent_y + attempts += 1 + # if attempts == attempts_tot: + # print("Could not place speckle without overlap, placing anyway.") + # Place the speckle finally + pixelsInDisk(cent_x, cent_y, + screen_size_width, screen_size_height, + speckle_size, foreground_colour, image, check_overlap=False) + else: + check_overlap = False + for i in range(int(total_speckles)): + cent_x = np.random.uniform(0, screen_size_width) + cent_y = np.random.uniform(0, screen_size_height) + pixelsInDisk(cent_x, cent_y, + screen_size_width, screen_size_height, + speckle_size, foreground_colour, image, check_overlap) + results[i, 0] = i + 1 + results[i, 1] = 1 + results[i, 2] = 2 + results[i, 3] = cent_x + results[i, 4] = cent_y + + # Apply Gaussian blur + image_blur = np.copy(image) + image_blur = ndimage.gaussian_filter(image_blur, sigma=sigma) + + return image_blur, results + + +def generate_speckles_random_disks_grid(screen_size_width: int, screen_size_height: int, + feature_size_width: float, feature_size_height: float, + foreground_colour: int, background_colour: int, + bit_depth: int, type_gen: str, seed: int, **kwargs) -> np.ndarray: + """A function to generate a speckle pattern image by randomly perturbating a grid of regularly-placed + disk-shaped speckles based on disrete uniform probability distribution. + + Parameters + ---------- + screen_size_width : int + Dimension of the image width-wise (pixels) + screen_size_height : int + Dimension of the image height-wise (pixels) + feature_size_width : float + Speckle size width-wise + feature_size_height : float + Speckle size height-wise + foreground_colour : int + Colour value for the speckles + background_colour : int + Colour value for the background + bit_depth : int + Bit depth of the image (8 or 16) + type_gen : str + Type of noise to generate (must be 'random_disks_grid' for this function) + seed : int + Random seed for the noise generation + total_speckles : int + Total number of speckles to place + sigma : float + Standard deviation for Gaussian blur applied after speckle placement + perturbation_max : float + Maximum amount to move speckles by during grid perturbation + + Returns + ------- + np.ndarray + Speckle pattern image as a 2D numpy array and speckle generation stats + """ + + speckle_size = (feature_size_width + feature_size_height) / 2 + total_speckles = kwargs.get("total_speckles") + sigma = kwargs.get("sigma") + perturbation_max = kwargs.get("perturbation_max") + + # check that none of the required parameters are None + assert speckle_size is not None, "speckle_size parameter is required." + assert total_speckles is not None, "total_speckles parameter is required." + assert sigma is not None, "sigma parameter is required." + assert type_gen == "random_disks_grid", "Type_gen must be 'random_disks_grid' for this function." + + np.random.seed(seed) + image: np.ndarray = np.ones((screen_size_height,screen_size_width), + dtype=np.uint16 if bit_depth == 16 else np.uint8) * background_colour + + results = np.zeros([total_speckles, 5]) # speckle number, attempts, with/without/not checked overlap (1/0/2), cent_x, cent_y + + # Calculate the size of each grid cell and create initial speckle positions + grid_cols = int(np.sqrt(total_speckles)) + grid_rows = int(np.ceil(total_speckles / grid_cols)) + grid_cell_width = screen_size_width // grid_cols + grid_cell_height = screen_size_height // grid_rows + + i = np.arange(grid_rows) + j = np.arange(grid_cols) + speckle_positions_x = j * grid_cell_width + grid_cell_width // 2 + speckle_size // 2 + speckle_positions_y = i * grid_cell_height + grid_cell_height // 2 + speckle_size // 2 + speckle_positions = np.array(np.meshgrid(speckle_positions_x, speckle_positions_y)).T.reshape(-1, 2) + + # Remove extra speckle positions randomly + total_speckles_diff = len(speckle_positions) - total_speckles + + if total_speckles_diff > 0: + remove_indices = np.random.choice(len(speckle_positions), total_speckles_diff, replace=False) + mask = np.ones(len(speckle_positions), dtype=bool) + mask[remove_indices] = False + speckle_positions = speckle_positions[mask] + + # Perturb each speckle's position + perturbations = np.random.randint(-perturbation_max, perturbation_max, size=(total_speckles, 2)) + perturbed_positions = [] + check_overlap = False + for idx, (x, y) in enumerate(speckle_positions): + # Add random perturbation to the x and y coordinates + x_offset, y_offset = perturbations[idx] + new_x = np.clip(x + x_offset, 0, screen_size_width - 1) + new_y = np.clip(y + y_offset, 0, screen_size_height - 1) + perturbed_positions.append((new_x, new_y)) + + pixelsInDisk(new_x, new_y, + screen_size_width, screen_size_height, + speckle_size, foreground_colour, image, check_overlap) + + results[idx, 0] = idx + 1 + results[idx, 1] = 1 + results[idx, 2] = 2 + results[idx, 3] = new_x + results[idx, 4] = new_y + + # Apply Gaussian blur + image_blur = np.copy(image) + image_blur = ndimage.gaussian_filter(image_blur, sigma=sigma) + + return image_blur, results + +def generate_speckles_perlin_noise(screen_size_width: int, screen_size_height: int, + feature_size_width: float, feature_size_height: float, + foreground_colour: int, background_colour: int, + bit_depth: int, type_gen: str, seed: int, **kwargs) -> np.ndarray: + """A function to generate a speckle pattern image using Perlin noise or fractal noise. + This function uses noise implementation from perlin-numpy package + (https://github.com/pvigier/perlin-numpy/tree/master) + + Parameters + ---------- + screen_size_width : int + Dimension of the image width-wise (pixels) + screen_size_height : int + Dimension of the image height-wise (pixels) + feature_size_width : float + Speckle size width-wise + feature_size_height : float + Speckle size height-wise + foreground_colour : int + Colour value for the speckles + background_colour : int + Colour value for the background + bit_depth : int + Bit depth of the image (8 or 16) + type_gen : str + Type of noise to generate (must be 'perlin' or 'fractal' for this function) + seed : int + Random seed for the noise generation + lacunarity : float (optional) + Lacunarity parameter for fractal noise (required if type_gen is 'fractal') + octaves : int (optional) + Number of octaves for fractal noise (required if type_gen is 'fractal') + + Returns + ------- + np.ndarray + Speckle pattern image as a 2D numpy array + """ + + res_width = int(screen_size_width / feature_size_width) + res_height = int(screen_size_height / feature_size_height) + + np.random.seed(seed) + + assert type_gen in ["perlin", "fractal"], "type_gen must be either 'perlin' or 'fractal' for this function." + + if type_gen == "perlin": + assert screen_size_width % res_width == 0, "The screen width must be a multiple of res_width. res_width = int(screen_size_width / feature_size_width)." + assert screen_size_height % res_height == 0, "The screen height must be a multiple of res_height. res_height = int(screen_size_height / feature_size_height)." + image = generate_perlin_noise_2d(shape = (screen_size_height, screen_size_width), + res = (res_height, res_width)) + + elif type_gen == "fractal": + lacunarity = kwargs.get("lacunarity", None) + octaves = kwargs.get("octaves", None) + assert lacunarity is not None, "lacunarity parameter is required for fractal noise." + assert octaves is not None, "octaves parameter is required for fractal noise." + assert screen_size_width % (lacunarity ** (octaves - 1) * res_width) == 0, "The screen width must be a multiple of lacunarity^(octaves-1)*res_width. res_width = int(screen_size_width / feature_size_width)." + assert screen_size_height % (lacunarity ** (octaves - 1) * res_height) == 0, "The screen height must be a multiple of lacunarity^(octaves-1)*res_height. res_height = int(screen_size_height / feature_size_height)." + image = generate_fractal_noise_2d(shape = (screen_size_height, screen_size_width), + res = (res_height, res_width), + octaves = octaves) + + # scale to background and foreground colours + min_val = np.min(image) + max_val = np.max(image) + image = (image - min_val) / (max_val - min_val) # Normalise to [0, 1] + + if foreground_colour >= background_colour: + # Black speckles on white background + image = image * (foreground_colour - background_colour) + background_colour + else: + # White speckles on black background + image = (1 - image) * (background_colour - foreground_colour) + foreground_colour + + image = image.astype(np.uint16 if bit_depth == 16 else np.uint8) + + return image + +def generate_speckles_simplex_noise(screen_size_width: int, screen_size_height: int, + feature_size_width: float, feature_size_height: float, + foreground_colour: int, background_colour: int, + bit_depth: int, type_gen: str, seed: int, **kwargs) -> np.ndarray: + """A function to generate a speckle pattern image using Simplex noise. + This function uses noise implementation from opensimplex package + (https://pypi.org/project/opensimplex/) + (https://code.larus.se/lmas/opensimplex) + + Parameters + ---------- + screen_size_width : int + Dimension of the image width-wise (pixels) + screen_size_height : int + Dimension of the image height-wise (pixels) + feature_size_width : float + Speckle size width-wise + feature_size_height : float + Speckle size height-wise + foreground_colour : int + Colour value for the speckles + background_colour : int + Colour value for the background + bit_depth : int + Bit depth of the image (8 or 16) + type_gen : str + Type of noise to generate (must be 'simplex' for this function) + seed : int + Random seed for the noise generation + + Returns + ------- + np.ndarray + Speckle pattern image as a 2D numpy array + """ + + assert type_gen == "simplex", "Type_gen must be 'simplex' for this function." + simplex.seed(seed) + + ix, iy = np.arange(screen_size_width), np.arange(screen_size_height) + image = simplex.noise2array(ix/feature_size_width, iy/feature_size_height) + + # scale to background and foreground colours + min_val = np.min(image) + max_val = np.max(image) + image = (image - min_val) / (max_val - min_val) # Normalise to [0, 1] + + if foreground_colour >= background_colour: + # Black speckles on white background + image = image * (foreground_colour - background_colour) + background_colour + else: + # White speckles on black background + image = (1 - image) * (background_colour - foreground_colour) + foreground_colour + + image = image.astype(np.uint16 if bit_depth == 16 else np.uint8) + + return image + +def generate_speckles(screen_size_width: int, screen_size_height: int, + feature_size_width: float, feature_size_height: float, + foreground_colour: int, background_colour: int, + bit_depth: int, type_gen: str, + seed: int, **kwargs) -> np.ndarray: + """A function to generate a speckle pattern image using specified generation method. + Speckle generation methods include 'random_disks', "random_disks_grid", 'perlin', 'fractal', and 'simplex'. + + Parameters + ---------- + screen_size_width : int + Dimension of the image width-wise (pixels) + screen_size_height : int + Dimension of the image height-wise (pixels) + feature_size_width : float + Speckle size width-wise + feature_size_height : float + Speckle size height-wise + foreground_colour : int + Colour value for the speckles + background_colour : int + Colour value for the background + bit_depth : int + Bit depth of the image (8 or 16) + type_gen : str + Type of noise to generate (must be one of 'random_disks', "random_disks_grid", 'perlin', 'fractal', 'simplex') + seed : int + Random seed for the noise generation + + Returns + ------- + np.ndarray + Speckle pattern image as a 2D numpy array and speckle generation stats (if applicable) + + Raises + ------ + ValueError + Unknown speckle generation type + """ + + assert bit_depth in [8, 16], "Bit depth should be either 8 or 16." + dispatch = { + "random_disks": generate_speckles_random_disks, + "random_disks_grid": generate_speckles_random_disks_grid, + "perlin": generate_speckles_perlin_noise, + "fractal": generate_speckles_perlin_noise, + "simplex": generate_speckles_simplex_noise + } + + generate = dispatch.get(type_gen) + if generate: + return generate(screen_size_width, screen_size_height, + feature_size_width, feature_size_height, + foreground_colour, background_colour, + bit_depth, type_gen, seed, **kwargs) + else: + raise ValueError(f"Unknown speckle generation type: {type_gen}") + diff --git a/src/pyvale/specklegen/speckle_pattern_diagnostics.py b/src/pyvale/specklegen/speckle_pattern_diagnostics.py new file mode 100644 index 00000000..58977444 --- /dev/null +++ b/src/pyvale/specklegen/speckle_pattern_diagnostics.py @@ -0,0 +1,301 @@ +import numpy as np +import matplotlib.pyplot as plt +import matplotlib +import seaborn as sns +import scipy.fftpack +from scipy.optimize import curve_fit +from scipy.stats import skew, kurtosis +from skimage.measure import shannon_entropy + + +def speckle_pattern_statistics(image: np.ndarray, dynamic_range: int) -> dict: + """A function to perform diagnostics on the speckle pattern image. + + Parameters + ---------- + image : np.ndarray + 2D numpy array representing the speckle pattern + dynamic_range : int + Maximum pixel value (irradiance) based on bit depth + + Returns + ------- + dict + Dictionary with diagnostic results + """ + + HFWHM, HeSquared, H_fit_stats, VFWHM, VeSquared, V_fit_stats, popt_H, popt_V, h_profile, v_profile = speckle_size(image) + avg_speckle_size_fwhm = np.mean([HFWHM, VFWHM]) + avg_speckle_size_e2 = np.mean([HeSquared, VeSquared]) + + # Calculate black/white ratio + black_pixels = np.sum(image < (dynamic_range // 2)) + white_pixels = np.sum(image >= (dynamic_range // 2)) + ratio = black_pixels / white_pixels + + # Calculate mean intensity gradient (simple finite difference) + grad_y, grad_x = np.gradient(image) + mean_gradient = np.mean(np.sqrt(grad_x ** 2 + grad_y ** 2)) + std_dev = np.std(image) + avg = np.mean(image) + contrast = std_dev / avg + skewness = skew(image.flatten()) + kurt = kurtosis(image.flatten()) + entropy = shannon_entropy(image) + peak_to_mean = np.max(image) / np.mean(image) + + results = { + "black_white_ratio": ratio, + "mean_intensity_gradient": mean_gradient, + "std_dev_irradiance": std_dev, + "avg_irradiance": avg, + "contrast": contrast, + "skewness": skewness, + "kurtosis": kurt, + "shannon_entropy": entropy, + "peak_to_mean_ratio": peak_to_mean, + "avg_speckle_size_fwhm": avg_speckle_size_fwhm, + "avg_speckle_size_e2": avg_speckle_size_e2, + "H_fit_stats": H_fit_stats, + "V_fit_stats": V_fit_stats + } + + return results + +def speckle_pattern_plots(image: np.ndarray, dynamic_range: int, + save_path: str = None) -> dict: + """A function to generate and save diagnostic plots for the speckle pattern. + + Parameters + ---------- + image : np.ndarray + 2D numpy array representing the speckle pattern + dynamic_range : int + Maximum pixel value (irradiance) based on bit depth + save_path : str (optional) + Path to save the generated plots using example formatting to the folder if provided + + Returns + ------- + dict + Dictionary containing figures and axes of the generated plots + """ + + HFWHM, HeSquared, H_fit_stats, VFWHM, VeSquared, V_fit_stats, popt_H, popt_V, h_profile, v_profile = speckle_size(image) + plots = {} + + # Set Seaborn theme + sns.set_theme(style="darkgrid") + matplotlib.rcParams['font.family'] = 'Sans-serif' + plt.rc('text', usetex=True) + fontsize1= 17 + + fig, axes = plt.subplots(1, 1, figsize=(7.0, 5.0)) + ax = axes + ax.imshow(image, cmap='gray', vmin=0, vmax=dynamic_range) + cbar = plt.colorbar(ax.images[0], ax=ax, fraction=0.046, pad=0.04) + cbar.ax.tick_params(labelsize=fontsize1-2) + cbar.set_label('Irradiance', fontsize=fontsize1-2) + ax.tick_params(axis='both', which='major', labelsize=fontsize1-2) + ax.set_xlabel("Position [pixel]", fontsize=fontsize1-2) + ax.set_ylabel("Position [pixel]", fontsize=fontsize1-2) + ax.set_title("Speckle pattern", fontsize=fontsize1) + + if save_path is not None: + plt.savefig(f"{save_path}/speckle_pattern.jpg", dpi=300, format='jpg', bbox_inches='tight') + + plots['speckle_pattern_fig'] = fig + plots['speckle_pattern_ax'] = ax + + f_image = scipy.fftpack.fft2(image) + f_image_shifted = scipy.fftpack.fftshift(f_image) + magnitude_spectrum = np.abs(f_image_shifted) + magnitude_spectrum_log = np.log1p(magnitude_spectrum) + fig, axes = plt.subplots(1, 1, figsize=(7.0, 5.0)) + ax = axes + ax.imshow(magnitude_spectrum_log, cmap='gray') + cbar = plt.colorbar(ax.images[0], ax=ax, fraction=0.046, pad=0.04) + cbar.ax.tick_params(labelsize=fontsize1-2) + ax.tick_params(axis='both', which='major', labelsize=fontsize1-2) + ax.set_title("Spatial frequency (log scale)", fontsize=fontsize1) + ax.set_xlabel("Frequency [1/pixel]", fontsize=fontsize1-2) + ax.set_ylabel("Frequency [1/pixel]", fontsize=fontsize1-2) + + if save_path is not None: + plt.savefig(f"{save_path}/frequency_spectrum.jpg", dpi=300, format='jpg', bbox_inches='tight') + + plots['frequency_spectrum_fig'] = fig + plots['frequency_spectrum_ax'] = ax + + fig, axes = plt.subplots(1, 1, figsize=(7.0, 5.0)) + ax = axes + ax.hist(image.ravel(), density=True, bins=int(dynamic_range/10), color='blue', alpha=0.7, log=True) + ax.set_title("Histogram of irradiance values", fontsize=fontsize1) + ax.set_xlabel("Pixel value", fontsize=fontsize1-2) + ax.set_ylabel("Density (log scale)", fontsize=fontsize1-2) + ax.tick_params(axis='both', which='major', labelsize=fontsize1-2) + + if save_path is not None: + plt.savefig(f"{save_path}/pixel_value_histogram.jpg", dpi=300, format='jpg', bbox_inches='tight') + + plots['pixel_value_histogram_fig'] = fig + plots['pixel_value_histogram_ax'] = ax + + plt.figure(figsize=(7.0, 5.0)) + x_H = np.arange(1, h_profile.size + 1) + x_V = np.arange(1, v_profile.size + 1) + plt.subplot(2, 1, 1) + plt.plot(x_H, h_profile, 'b-', label='Horisontal autocov.', + linewidth=2) + plt.plot(x_H, gaussian(x_H, *popt_H), 'r--', label='Gaussian interpol.', linewidth=2) + plt.title('Horisontal autocovariance', fontsize=fontsize1) + plt.xlabel('Lag [pixels]', fontsize=fontsize1-2) + plt.ylabel('Autocov. ' + r"[pixel$^2$]", fontsize=fontsize1-2) + plt.legend(fontsize=fontsize1-2) + plt.tick_params(axis='both', which='major', labelsize=fontsize1-2) + plt.subplot(2, 1, 2) + plt.plot(x_V, v_profile, 'b-', label='Vertical autocov.', + linewidth=2) + plt.plot(x_V, gaussian(x_V, *popt_V), 'r--', label='Gaussian interpol.', linewidth=2) + plt.title('Vertical autocovariance', fontsize=fontsize1) + plt.xlabel('Lag [pixels]', fontsize=fontsize1- 2) + plt.ylabel('Autocov. ' + r"[pixel$^2$]", fontsize=fontsize1-2) + plt.legend(fontsize=fontsize1-2) + plt.tick_params(axis='both', which='major', labelsize=fontsize1-2) + plt.tight_layout() + + if save_path is not None: + plt.savefig(f"{save_path}/autocovariance.jpg", dpi=300, format='jpg', bbox_inches='tight') + + # Extract figure and axis + fig = plt.gcf() + ax = plt.gca() + + plots['autocovariance_fig'] = fig + plots['autocovariance_ax'] = ax + + return plots + +def speckle_size(image: np.ndarray) -> tuple: + """ A function to calculate speckle size from the autocovariance of the speckle pattern. + + Parameters + ---------- + image : np.ndarray + 2D numpy array representing the speckle pattern + + Returns + ------- + tuple + HFWHM : float + Full width at half maximum for horisontal profile. + HeSquared : float + 1/e^2 width for horisontal profile. + H_fit_stats : dict + R-squared goodness of fit for horisontal profile. + VFWHM : float + Full width at half maximum for vertical profile. + VeSquared : float + 1/e^2 width for vertical profile. + V_fit_stats : dict + R-squared goodness of fit for vertical profile. + popt_H : array + Optimal parameters for horisontal Gaussian fit. + popt_V : array + Optimal parameters for vertical Gaussian fit. + h_profile : np.ndarray + Horisontal autocovariance profile. + v_profile : np.ndarray + Vertical autocovariance profile. + """ + # image = (image - np.mean(image)) / np.std(image) + image = (image - np.mean(image)) + + # Autocorrelation (FFT analysis) + f_image_norm = scipy.fftpack.fft2(image) + power_spectrum = np.abs(f_image_norm)**2 + autocorr = scipy.fftpack.ifft2(power_spectrum).real + autocorr = scipy.fftpack.fftshift(autocorr) + autocorr /= np.max(autocorr) + centre_y, centre_x = np.array(autocorr.shape) // 2 + h_profile = autocorr[centre_y, :] + v_profile = autocorr[:, centre_x] + + HFWHM, HeSquared, H_fit_stats, VFWHM, VeSquared, V_fit_stats, popt_H, popt_V = fit_gaussian(h_profile, v_profile) + + return HFWHM, HeSquared, H_fit_stats, VFWHM, VeSquared, V_fit_stats, popt_H, popt_V, h_profile, v_profile + +def gaussian(x: np.ndarray, a1: float, b1: float, c1: float) -> np.ndarray: + """A function to model Gaussian distribution. + + Parameters + ---------- + x : np.ndarray + Input array. + a1 : float + Peak amplitude (height) + b1 : float + Center position (mean) + c1 : float + Standard deviation (width, spread) + + Returns + ------- + np.ndarray + _description_ + """ + return a1 * np.exp(-((x - b1)/c1) ** 2) + +def fit_gaussian(H: np.ndarray, V: np.ndarray) -> tuple: + """ Fit Gaussian functions to the horisontal and vertical autocovariance profiles. + + Parameters + ---------- + H : np.ndarray + Horisontal autocovariance profile. + V : np.ndarray + Vertical autocovariance profile. + + Returns + ------- + tuple + HFWHM : float + Full width at half maximum for horisontal profile. + HeSquared : float + 1/e^2 width for horisontal profile. + H_fit_stats : dict + R-squared goodness of fit for horisontal profile. + VFWHM : float + Full width at half maximum for vertical profile. + VeSquared : float + 1/e^2 width for vertical profile. + V_fit_stats : dict + R-squared goodness of fit for vertical profile. + popt_H : array + Optimal parameters for horisontal Gaussian fit. + popt_V : array + Optimal parameters for vertical Gaussian fit. + """ + + range_H = np.arange(1, H.size + 1) + range_V = np.arange(1, V.size + 1) + + low_H = np.where(H > 0.2)[0] + low_V = np.where(V > 0.2)[0] + + popt_H, _ = curve_fit(gaussian, range_H[low_H], H[low_H], p0=[1, np.argmax(H), 1]) + popt_V, _ = curve_fit(gaussian, range_V[low_V], V[low_V], p0=[1, np.argmax(V), 1]) + + # Extract FWHM and 1/e^2 widths + HFWHM = 2 * popt_H[2] * np.sqrt(-np.log(0.5 / popt_H[0])) # FWHM for horisontal + VFWHM = 2 * popt_V[2] * np.sqrt(-np.log(0.5 / popt_V[0])) # FWHM for vertical + HeSquared = popt_H[2] * np.sqrt(-np.log(0.1353353 / popt_H[0])) # 1/e^2 for horisontal + VeSquared = popt_V[2] * np.sqrt(-np.log(0.1353353 / popt_V[0])) # 1/e^2 for vertical + + # R-squared goodness of fit + H_fit_stats = {'R_squared': 1 - np.sum((H[low_H] - gaussian(range_H[low_H], *popt_H)) ** 2) / np.sum((H[low_H] - np.mean(H[low_H])) ** 2)} + V_fit_stats = {'R_squared': 1 - np.sum((V[low_V] - gaussian(range_V[low_V], *popt_V)) ** 2) / np.sum((V[low_V] - np.mean(V[low_V])) ** 2)} + + return HFWHM, HeSquared, H_fit_stats, VFWHM, VeSquared, V_fit_stats, popt_H, popt_V + + diff --git a/src/pyvale/verif/specklegenconst.py b/src/pyvale/verif/specklegenconst.py new file mode 100644 index 00000000..6c35f21a --- /dev/null +++ b/src/pyvale/verif/specklegenconst.py @@ -0,0 +1,18 @@ +#=============================================================================== +# pyvale: the python validation engine +# License: MIT +# Copyright (C) 2025 The Computer Aided Validation Team +#=============================================================================== +from pathlib import Path + +""" +DEVELOPER VERIFICATION MODULE +-------------------------------------------------------------------------------- +This module contains developer utility functions used for verification testing +of the speckle generation pattern in pyvale. + +Specifically, this module contains constants used for verification testing. +""" + +GOLD_PATH: Path = Path.cwd() / "tests" / "specklegen" / "gold" +GOLD_SEED: int = 123 diff --git a/src/pyvale/verif/specklegold.py b/src/pyvale/verif/specklegold.py new file mode 100644 index 00000000..b73e53ec --- /dev/null +++ b/src/pyvale/verif/specklegold.py @@ -0,0 +1,99 @@ +#=============================================================================== +# pyvale: the python validation engine +# License: MIT +# Copyright (C) 2025 The Computer Aided Validation Team +#=============================================================================== + +""" +DEVELOPER VERIFICATION MODULE +-------------------------------------------------------------------------------- +This module contains developer utility functions used for verification testing +of the speckle generation pattern in pyvale. + +Specifically, this module contains gold measurement generation function. +""" + +import numpy as np +import os +import pyvale.verif.specklegenconst as specklegenconst +import pyvale.specklegen as specklegen + + +def gen_gold_measurements(param_dict: dict) -> None: + case_tot = param_dict.get("case_tot") + type_gen = param_dict.get("type_gen") + + screen_size_width = param_dict.get("screen_size_width") + screen_size_height = param_dict.get("screen_size_height") + speckle_size = param_dict.get("speckle_size") + bit_depth = param_dict.get("bit_depth") + theme = param_dict.get("theme") + seed = param_dict.get("seed") + reduce_overlap = param_dict.get("reduce_overlap") + sigma = param_dict.get("sigma") + attempts_tot = param_dict.get("attempts_tot") + perturbation_max = param_dict.get("perturbation_max") + octaves = param_dict.get("octaves") + lacunarity = param_dict.get("lacunarity") + feature_size_width = speckle_size + feature_size_height = speckle_size + + speckle_area = np.pi * (speckle_size / 2) ** 2 + total_area = screen_size_width * screen_size_height + total_speckles = int((0.5 * total_area) / speckle_area) + dynamic_range: int = 2**bit_depth - 1 + background_colour = 0 if theme == 'white_on_black' else dynamic_range + foreground_colour = dynamic_range if theme == 'white_on_black' else 0 + + + for i in range(case_tot): + print(f"Generating gold output for case: {type_gen}_{i+1}") + + if type_gen == "random_disks": + image, _ = specklegen.generate_speckles(screen_size_width, screen_size_height, + feature_size_width, feature_size_height, + foreground_colour, background_colour, + bit_depth, type_gen, seed, + total_speckles=total_speckles, + reduce_overlap=reduce_overlap, + sigma=sigma, attempts_tot=attempts_tot, + perturbation_max=perturbation_max) + elif type_gen == "random_disks_grid": + image, _ = specklegen.generate_speckles(screen_size_width, screen_size_height, + feature_size_width, feature_size_height, + foreground_colour, background_colour, + bit_depth, type_gen, seed, + total_speckles=total_speckles, + reduce_overlap=reduce_overlap, + sigma=sigma, perturbation_max=perturbation_max) + elif type_gen == "perlin": + image = specklegen.generate_speckles(screen_size_width, screen_size_height, + feature_size_width, feature_size_height, + foreground_colour, background_colour, + bit_depth, type_gen, seed) + elif type_gen == "fractal": + image = specklegen.generate_speckles(screen_size_width, screen_size_height, + feature_size_width, feature_size_height, + foreground_colour, background_colour, + bit_depth, type_gen, seed, + octaves=octaves, lacunarity=lacunarity) + + elif type_gen == "simplex": + image = specklegen.generate_speckles(screen_size_width, screen_size_height, + feature_size_width, feature_size_height, + foreground_colour, background_colour, + bit_depth, type_gen, seed) + + save_path = specklegenconst.GOLD_PATH + if not os.path.exists(save_path): + os.makedirs(save_path) + np.save(f"{save_path}/{type_gen}_image_{i+1}.npy", image) + + if i == 0: + image_ref = image + else: + error_abs = np.abs(image - image_ref) + np.save(f"{save_path}/{type_gen}_image_abs_error_{i+1}.npy", image) + print(f"Avg. absolute error {np.mean(error_abs)}, max. absolute error {np.max(error_abs)}") + + diff --git a/tests/specklegen/generate_speckles_test.py b/tests/specklegen/generate_speckles_test.py new file mode 100644 index 00000000..5b9cf1e6 --- /dev/null +++ b/tests/specklegen/generate_speckles_test.py @@ -0,0 +1,183 @@ +import numpy as np +import argparse +import json +import pytest +import sys +import os +sys.path.append(os.path.abspath(os.path.join(os.path.dirname(__file__), '..'))) +import pyvale.specklegen as specklegen + +width_range = np.arange(10, 1000, 250) +height_range = np.arange(10, 1000, 250) +speckle_size_range = np.arange(4, 25, 5) +total_speckles_range = np.arange(1, 10, 1) +sizes = list(zip(width_range, height_range, speckle_size_range)) + +print(width_range.shape) +print(height_range.shape) +print(speckle_size_range.shape) + +pytestmark = pytest.mark.parametrize("width,height,speckle_size", sizes) + +def test_pixelsInDisk_out_of_bounds(width: int, height: int, speckle_size: int) -> None: + img = np.zeros((height, width), dtype=np.uint8) + cantre_x, centre_y = -speckle_size * 2, -speckle_size * 2 # Out of bounds + fg_colour = 1 + check_overlap = False + + specklegen.pixelsInDisk(cantre_x, centre_y, width, height, speckle_size, fg_colour, img, check_overlap) + + assert np.all(img == 0) + +def test_pixelsInDisk_basic_circle(width: int, height: int, speckle_size: int) -> None: + img = np.zeros((height, width), dtype=np.uint8) + cantre_x, centre_y = 5, 5 + fg_colour = 1 + check_overlap = False + + specklegen.pixelsInDisk(cantre_x, centre_y, width, height, speckle_size, fg_colour, img, check_overlap) + + # Check that the centre pixel is set + assert img[centre_y, cantre_x] == fg_colour + + # Check that pixels outside the radius are not set + radius_sq = (speckle_size / 2) ** 2 + for y in range(height): + for x in range(width): + dist_sq = (x + 0.5 - cantre_x) ** 2 + (y + 0.5 - centre_y) ** 2 + if dist_sq > radius_sq: + assert img[y, x] == 0 + +def test_pixelsInDisk_near_edges(width: int, height: int, speckle_size: int) -> None: + img = np.zeros((height, width), dtype=np.uint8) + cantre_x, centre_y = 0, 0 + fg_colour = 1 + check_overlap = False + + specklegen.pixelsInDisk(cantre_x, centre_y, width, height, speckle_size, fg_colour, img, check_overlap) + + # Pixels outside bounds should not be accessed => Mo error => Passed test + assert img[0, 0] == fg_colour # Centre pixel should be set + + radius_sq = (speckle_size / 2) ** 2 + for y in range(height): + for x in range(width): + dist_sq = (x + 0.5 - cantre_x) ** 2 + (y + 0.5 - centre_y) ** 2 + if dist_sq <= radius_sq: + assert img[y, x] == fg_colour + +@pytest.mark.parametrize("total_speckles", total_speckles_range) +def test_generate_speckles_basic_no_overlap(width: int, height: int, speckle_size: int, + total_speckles: int) -> None: + fg_colour = 255 + bg_colour = 0 + total_speckles = 5 + sigma = 1.0 + bit_depth = 8 + feature_size_width = speckle_size + feature_size_height = speckle_size + seed = 123 + type_gens = ["random_disks", "random_disks_grid"] + reduce_overlap=False + perturbation_max = 6 + + for type_gen in type_gens: + image, results = specklegen.generate_speckles(width, height, + feature_size_width, feature_size_height, + fg_colour, bg_colour, + bit_depth, type_gen, seed, + total_speckles=total_speckles, + reduce_overlap=reduce_overlap, + sigma=sigma, perturbation_max=perturbation_max) + + assert image.shape == (height, width) + assert results.shape == (total_speckles, 5) + + # All speckles should have attempt count = 1 and overlap flag = 2 (not checked) + assert np.all(results[:, 1] == 1) + assert np.all(results[:, 2] == 2) + + +@pytest.mark.parametrize("total_speckles", total_speckles_range) +def test_generate_speckles_with_overlap_reduction(width: int, height: int, speckle_size: int, + total_speckles: int) -> None: + fg_colour = 255 + bg_colour = 0 + sigma = 1.0 + bit_depth = 8 + feature_size_width = speckle_size + feature_size_height = speckle_size + seed = 123 + type_gen = "random_disks" + reduce_overlap=True + attempts_tot=50 + + image, results = specklegen.generate_speckles(width, height, + feature_size_width, feature_size_height, + fg_colour, bg_colour, + bit_depth, type_gen, seed, + total_speckles=total_speckles, + reduce_overlap=reduce_overlap, + sigma=sigma, attempts_tot=attempts_tot) + + assert image.shape == (height, width) + assert results.shape == (total_speckles, 5) + + # Overlap flag should be either 0 or 1 (placed with or without overlap) + assert np.all(np.isin(results[:, 2], [0, 1])) + +@pytest.mark.parametrize("total_speckles", total_speckles_range) +def test_generate_speckles_bit_depth(width: int, height: int, speckle_size: int, + total_speckles: int) -> None: + fg_colour = 65535 + bg_colour = 0 + sigma = 0 + bit_depth = 16 + feature_size_width = speckle_size + feature_size_height = speckle_size + seed = 123 + type_gens = ["random_disks", "random_disks_grid"] + reduce_overlap=False + perturbation_max = 6 + + for type_gen in type_gens: + image, results = specklegen.generate_speckles(width, height, + feature_size_width, feature_size_height, + fg_colour, bg_colour, + bit_depth, type_gen, seed, + total_speckles=total_speckles, + reduce_overlap=reduce_overlap, + sigma=sigma, perturbation_max=perturbation_max) + + assert image.dtype == np.uint16 + assert image.max() <= fg_colour + assert results.shape == (total_speckles, 5) + +@pytest.mark.parametrize("total_speckles", total_speckles_range) +def test_generate_speckles_blur(width: int, height: int, speckle_size: int, + total_speckles: int) -> None: + fg_colour = 255 + bg_colour = 0 + sigma = 1.5 + bit_depth = 8 + feature_size_width = speckle_size + feature_size_height = speckle_size + seed = 123 + type_gens = ["random_disks", "random_disks_grid"] + reduce_overlap=False + perturbation_max = 6 + + for type_gen in type_gens: + image, results = specklegen.generate_speckles(width, height, + feature_size_width, feature_size_height, + fg_colour, bg_colour, + bit_depth, type_gen, seed, + total_speckles=total_speckles, + reduce_overlap=reduce_overlap, + sigma=sigma, perturbation_max=perturbation_max) + + assert np.any(image != 0) + assert np.any(image != fg_colour) + assert image.shape == (height, width) + assert results.shape == (total_speckles, 5) +