From c957ca84e60a7a0cbb63ade96c386a7b00a4aa31 Mon Sep 17 00:00:00 2001 From: Jamy Date: Tue, 21 Mar 2023 15:19:14 +0100 Subject: [PATCH 1/4] Create fast_monte_carlo.py --- fast_monte_carlo.py | 125 ++++++++++++++++++++++++++++++++++++++++++++ 1 file changed, 125 insertions(+) create mode 100644 fast_monte_carlo.py diff --git a/fast_monte_carlo.py b/fast_monte_carlo.py new file mode 100644 index 0000000..5e1b2f4 --- /dev/null +++ b/fast_monte_carlo.py @@ -0,0 +1,125 @@ +# -*- coding: utf-8 -*- +""" +Created on Tue Mar 21 11:49:06 2023 + +@author: jamyl +""" + +import numpy as np +from tqdm import tqdm +import os +import matplotlib.pyplot as plt + +n_patches = int(1e5) +n_brightness_levels = 1000 + +alpha = 1.80710882e-4 # from measurements +beta = 3.1937599182128e-6 # from https://www.photonstophotos.net/Charts/RN_ADU.htm chart +tol = 3 + +#%% + +def get_non_linearity_bound(alpha, beta, tol): + tol_sq = tol*tol + xmin = tol_sq/2 * (alpha + np.sqrt(tol_sq*alpha*alpha + 4*beta)) + + xmax = (2 + tol_sq*alpha - np.sqrt((2+tol_sq*alpha)**2 - 4*(1+tol_sq*beta)))/2 + + return xmin, xmax + + + +def unitary_MC(b): + # create the patch + patch = np.ones((n_patches, 3, 3)) * b + + # add noise and clip + patch1 = patch + np.sqrt(patch * alpha + beta) * np.random.randn(*patch.shape) + patch1 = np.clip(patch1, 0.0, 1.0) + + patch2 = patch + np.sqrt(patch * alpha + beta) * np.random.randn(*patch.shape) + patch2 = np.clip(patch2, 0.0, 1.0) + + # compute statistics and store + std_mean = 0.5 * np.mean((np.std(patch1, axis=(1, 2)) + np.std(patch2, axis=(1, 2)))) + + curr_mean1 = np.mean(patch1, axis=(1, 2)) + curr_mean2 = np.mean(patch2, axis=(1, 2)) + diff_mean = np.mean(np.abs(curr_mean1 - curr_mean2)) + + return diff_mean, std_mean + +def regular_MC(b_array): + sigmas = np.empty_like(b_array) + diffs = np.empty_like(b_array) + for i in tqdm(range(b_array.size)): + sigmas[i], diffs[i] = unitary_MC(b_array[i]) + return sigmas, diffs + +def interp_MC(b_array, + sigma_min, sigma_max, + diff_min, diff_max): + norm_b = (b_array - b_array[0])/(b_array[-1] - b_array[0]) + + sigmas_sq_lin = norm_b * (sigma_max**2 - sigma_min**2) + sigma_min**2 + diffs_sq_lin = norm_b * (diff_max**2 - diff_min**2) + diff_min**2 + + return np.sqrt(sigmas_sq_lin[1:-1]), np.sqrt(diffs_sq_lin[1:-1]) + + +def run_fast_MC(): + xmin, xmax = get_non_linearity_bound(alpha, beta, tol) + + # these 2 inddexes define the points for which the linear model is valid. + # The regular Monte carlo is applied to these points, that are then used as reference + # for the linear model + imin = int(np.ceil(xmin * n_brightness_levels)) + 1 + imax = int(np.floor(xmax * n_brightness_levels)) - 1 + + + sigmas = np.empty(n_brightness_levels+1) + diffs = np.empty(n_brightness_levels+1) + brigntess = np.arange(n_brightness_levels+1)/n_brightness_levels + + # running regular MC for non linear parts + nl_brigntess = np.concatenate((brigntess[:imin+1], brigntess[imax:])) + sigma_nl, diffs_nl = regular_MC(nl_brigntess) + sigmas[:imin+1], diffs[:imin+1] = sigma_nl[:imin+1], diffs_nl[:imin+1] + sigmas[imax:], diffs[imax:] = sigma_nl[imin+1:], diffs_nl[imin+1:] + + # padding using linear interpolation + brightness_l = brigntess[imin-1:imax+2] + + sigmas_l ,diffs_l = interp_MC(brightness_l, + sigmas[imin], sigmas[imax], + diffs[imin], diffs[imax]) + sigmas[imin:imax+1] = sigmas_l + diffs[imin:imax+1] = diffs_l + + + + + + plt.figure("sigma") + plt.plot(np.linspace(0, 1, n_brightness_levels+1), sigmas**2) + plt.grid() + plt.axline((xmin, 0), (xmin, 1e-5), c='k') + plt.axline((xmax, 0), (xmax, 1e-5), c='k') + plt.ylabel('Sigma**2') + plt.xlabel('brightness') + + plt.figure("Diff**2") + plt.plot(np.linspace(0, 1, n_brightness_levels+1), diffs**2) + plt.grid() + plt.axline((xmin, 0), (xmin, 1e-5), c='k') + plt.axline((xmax, 0), (xmax, 1e-5), c='k') + plt.ylabel('Diff**2') + plt.xlabel('brightness') + + + +run_fast_MC() + + + + \ No newline at end of file From fa9eabdbd422e03c3c4c90843c4352c09aaeaab2 Mon Sep 17 00:00:00 2001 From: Jamy Date: Fri, 24 Mar 2023 19:03:21 +0100 Subject: [PATCH 2/4] Added fast MC + .dng output --- fast_monte_carlo.py | 55 ++++--- handheld_super_resolution/super_resolution.py | 7 + handheld_super_resolution/utils_dng.py | 150 ++++++++++++++++++ 3 files changed, 194 insertions(+), 18 deletions(-) create mode 100644 handheld_super_resolution/utils_dng.py diff --git a/fast_monte_carlo.py b/fast_monte_carlo.py index 5e1b2f4..795bbbe 100644 --- a/fast_monte_carlo.py +++ b/fast_monte_carlo.py @@ -9,13 +9,18 @@ from tqdm import tqdm import os import matplotlib.pyplot as plt +import multiprocessing +from functools import partial +import time n_patches = int(1e5) n_brightness_levels = 1000 -alpha = 1.80710882e-4 # from measurements -beta = 3.1937599182128e-6 # from https://www.photonstophotos.net/Charts/RN_ADU.htm chart +# alpha = 1.80710882e-4 +# beta = 3.1937599182128e-6 tol = 3 +# While I +- tol*sigma is not saturated, the linear model is considered to hold. +# 3 and 5 are good values #%% @@ -29,7 +34,7 @@ def get_non_linearity_bound(alpha, beta, tol): -def unitary_MC(b): +def unitary_MC(alpha, beta, b): # create the patch patch = np.ones((n_patches, 3, 3)) * b @@ -49,11 +54,23 @@ def unitary_MC(b): return diff_mean, std_mean -def regular_MC(b_array): +def regular_MC(b_array, alpha, beta): + if multiprocessing.cpu_count() > 1: + N_CPU = multiprocessing.cpu_count()-1 + else: + N_CPU = 1 + multiprocessing.freeze_support() + pool = multiprocessing.Pool(processes=N_CPU) + func = partial(unitary_MC, alpha, beta) + + sigmas = np.empty_like(b_array) diffs = np.empty_like(b_array) - for i in tqdm(range(b_array.size)): - sigmas[i], diffs[i] = unitary_MC(b_array[i]) + for b_, result in enumerate(tqdm(pool.imap(func, list(b_array)), total=b_array.size,desc="Brightnesses")): + diffs[b_] = result[0] + sigmas[b_] = result[1] + + return sigmas, diffs def interp_MC(b_array, @@ -67,7 +84,9 @@ def interp_MC(b_array, return np.sqrt(sigmas_sq_lin[1:-1]), np.sqrt(diffs_sq_lin[1:-1]) -def run_fast_MC(): +def run_fast_MC(alpha, beta): + t0 = time.perf_counter() + print("Beginning...") xmin, xmax = get_non_linearity_bound(alpha, beta, tol) # these 2 inddexes define the points for which the linear model is valid. @@ -83,7 +102,7 @@ def run_fast_MC(): # running regular MC for non linear parts nl_brigntess = np.concatenate((brigntess[:imin+1], brigntess[imax:])) - sigma_nl, diffs_nl = regular_MC(nl_brigntess) + sigma_nl, diffs_nl = regular_MC(nl_brigntess, alpha, beta) sigmas[:imin+1], diffs[:imin+1] = sigma_nl[:imin+1], diffs_nl[:imin+1] sigmas[imax:], diffs[imax:] = sigma_nl[imin+1:], diffs_nl[imin+1:] @@ -95,31 +114,31 @@ def run_fast_MC(): diffs[imin], diffs[imax]) sigmas[imin:imax+1] = sigmas_l diffs[imin:imax+1] = diffs_l - + print('Finished! Estimated 1 noise curve in {:2f} sec'.format(time.perf_counter() - t0)) plt.figure("sigma") - plt.plot(np.linspace(0, 1, n_brightness_levels+1), sigmas**2) + plt.plot(np.linspace(0, 1, n_brightness_levels+1), sigmas) plt.grid() plt.axline((xmin, 0), (xmin, 1e-5), c='k') plt.axline((xmax, 0), (xmax, 1e-5), c='k') - plt.ylabel('Sigma**2') + plt.ylabel('Sigma') plt.xlabel('brightness') - plt.figure("Diff**2") - plt.plot(np.linspace(0, 1, n_brightness_levels+1), diffs**2) + plt.figure("Diff") + plt.plot(np.linspace(0, 1, n_brightness_levels+1), diffs) plt.grid() plt.axline((xmin, 0), (xmin, 1e-5), c='k') plt.axline((xmax, 0), (xmax, 1e-5), c='k') - plt.ylabel('Diff**2') + plt.ylabel('Diff') plt.xlabel('brightness') +if __name__ == "__main__": + run_fast_MC(alpha=1.80710882e-4, + beta=3.1937599182128e-6) -run_fast_MC() - - - + \ No newline at end of file diff --git a/handheld_super_resolution/super_resolution.py b/handheld_super_resolution/super_resolution.py index 6f847fc..4204c72 100644 --- a/handheld_super_resolution/super_resolution.py +++ b/handheld_super_resolution/super_resolution.py @@ -27,6 +27,7 @@ from . import raw2rgb from .utils import getTime, DEFAULT_NUMPY_FLOAT_TYPE, divide, add, round_iso from .utils_image import compute_grey_images, frame_count_denoising_gauss, frame_count_denoising_median +from .utils_dng import save_as_dng from .merge import merge, merge_ref from .kernels import estimate_kernels from .block_matching import init_block_matching, align_image_block_matching @@ -521,6 +522,12 @@ def process(burst_path, options=None, custom_params=None): #___ Running the handheld pipeline handheld_output, debug_dict = main(ref_raw.astype(DEFAULT_NUMPY_FLOAT_TYPE), raw_comp.astype(DEFAULT_NUMPY_FLOAT_TYPE), options, params) + #___ Saving output as dng if required + if 'dng_outpath' in options.keys(): + outpath = Path(options["dng_outpath"]) + if verbose_1: + print('Saving output to {}'.format(outpath.with_suffix('.dng').as_posix())) + save_as_dng(handheld_output.copy_to_host(), raw_path_list[ref_id], outpath) #___ Performing frame count aware denoising if enabled diff --git a/handheld_super_resolution/utils_dng.py b/handheld_super_resolution/utils_dng.py new file mode 100644 index 0000000..a076b6c --- /dev/null +++ b/handheld_super_resolution/utils_dng.py @@ -0,0 +1,150 @@ +# -*- coding: utf-8 -*- +""" +Created on Tue Mar 21 16:44:36 2023 + +@author: jamyl +""" + +import os + +import numpy as np +import exifread +import rawpy +import imageio + +EXIFTOOL_PATH = 'C:/Users/jamyl/Downloads/exiftool-12.44/exiftool.exe' +DNG_VALIDATE_PATH = 'C:/Users/jamyl/Downloads/dng_sdk_1_6/dng_sdk_1_6/dng_sdk/targets/win/debug64_x64/dng_validate.exe' + +def save_as_dng(np_img, ref_dng_path, outpath): + ''' + Saves a RGB numpy image as dng. + The image is first saved as 16bits tiff, then the extension is swapped + to .dng. The metadata are then overwritten using a reference dng, and + the final dng is built using dng_validate. + + Requires : + - dng_validate (can be found in dng sdk): + https://helpx.adobe.com/camera-raw/digital-negative.html#dng_sdk_download + + - exiftool + https://exiftool.org/ + + + Based on : + https://github.com/gluijk/dng-from-tiff/blob/main/dngmaker.bat + https://github.com/antonwolf/dng_stacker/blob/master/dng_stacker.bat + + Parameters + ---------- + np_img : numpy array + RGB image + rawpy_ref : _rawpy.Rawpy + image containing some relevant tags + outpath : Path + output save path. + + Returns + ------- + None. + + ''' + assert np_img.shape[-1] == 3 + + + np_int_img = np.copy(np_img) # copying to avoid inplace-overwritting + # Undo White balance and black level + ## get tags + with open(ref_dng_path, 'rb') as raw_file: + tags = exifread.process_file(raw_file) + + white_level = tags['Image Tag 0xC61D'].values[0] # there is only one white level + + black_levels = tags['Image BlackLevel'] + if isinstance(black_levels.values[0], int): + black_levels = np.array(black_levels.values) + else: # Sometimes this tag is a fraction object for some reason. It seems that black levels are all integers anyway + black_levels = np.array([int(x.decimal()) for x in black_levels.values]) + + raw = rawpy.imread(ref_dng_path) + white_balance = raw.camera_whitebalance + + ## Reverse WB + + for c in range(3): + np_int_img[:, :, c] /= white_balance[c] / white_balance[1] + np_int_img[:, :, c] = np_int_img[:, :, c] * (white_level - black_levels[c]) + black_levels[c] + + np_int_img = np.clip(np_int_img, 0, 2**16-1).astype(np.uint16) + + # Saving the image as 16 bits RGB tiff + save_as_tiff(np_int_img, outpath) + + tmp_path = outpath.parent/'tmp.dng' + + # Deleting tmp.dng if it is already existing + if os.path.exists(tmp_path): + os.remove(tmp_path) + + + + # Overwritting the tiff tags with dng tags, and replacing the .tif extension + # by .dng + cmd = ''' + {} -n\ + -IFD0:SubfileType#=0\ + -IFD0:PhotometricInterpretation#=34892\ + -SamplesPerPixel#=3\ + -overwrite_original -tagsfromfile {}\ + "-all:all>all:all"\ + -DNGVersion\ + -DNGBackwardVersion\ + -ColorMatrix1 -ColorMatrix2\ + "-IFD0:BlackLevelRepeatDim Date: Mon, 27 Mar 2023 14:09:00 +0200 Subject: [PATCH 3/4] Cleaned repo --- .../fast_monte_carlo.py | 141 ++++++-- handheld_super_resolution/super_resolution.py | 193 +++------- handheld_super_resolution/utils_dng.py | 192 ++++++++-- run_handheld.py | 337 +++++++++--------- 4 files changed, 503 insertions(+), 360 deletions(-) rename fast_monte_carlo.py => handheld_super_resolution/fast_monte_carlo.py (53%) diff --git a/fast_monte_carlo.py b/handheld_super_resolution/fast_monte_carlo.py similarity index 53% rename from fast_monte_carlo.py rename to handheld_super_resolution/fast_monte_carlo.py index 795bbbe..6b2aa98 100644 --- a/fast_monte_carlo.py +++ b/handheld_super_resolution/fast_monte_carlo.py @@ -2,12 +2,16 @@ """ Created on Tue Mar 21 11:49:06 2023 +This script is an extension of monte_carlo_simulation, optimised to run faster. +The simulation uses multiple CPU cores and focuses on the non linear parts of +the curve. Linear parts are interpolated, resulting in a considerable time gain. + + @author: jamyl """ import numpy as np from tqdm import tqdm -import os import matplotlib.pyplot as plt import multiprocessing from functools import partial @@ -18,8 +22,10 @@ # alpha = 1.80710882e-4 # beta = 3.1937599182128e-6 -tol = 3 -# While I +- tol*sigma is not saturated, the linear model is considered to hold. +TOL = 3 + +# While I +- tol*sigma is between 0 and 1, the linear model is considered +# to hold, and values are interpolated. # 3 and 5 are good values #%% @@ -35,6 +41,27 @@ def get_non_linearity_bound(alpha, beta, tol): def unitary_MC(alpha, beta, b): + """ + Runs a MC scheme to estimate sigma and d for a given brightness, alpha and + beta. + + Parameters + ---------- + alpha : TYPE + DESCRIPTION. + beta : TYPE + DESCRIPTION. + b : float in 0, 1 + brighntess + + Returns + ------- + diff_mean : float + mean difference + std_mean : float + mean standard deviation + + """ # create the patch patch = np.ones((n_patches, 3, 3)) * b @@ -45,7 +72,7 @@ def unitary_MC(alpha, beta, b): patch2 = patch + np.sqrt(patch * alpha + beta) * np.random.randn(*patch.shape) patch2 = np.clip(patch2, 0.0, 1.0) - # compute statistics and store + # compute statistics std_mean = 0.5 * np.mean((np.std(patch1, axis=(1, 2)) + np.std(patch2, axis=(1, 2)))) curr_mean1 = np.mean(patch1, axis=(1, 2)) @@ -55,6 +82,26 @@ def unitary_MC(alpha, beta, b): return diff_mean, std_mean def regular_MC(b_array, alpha, beta): + """ + Runs MC on multiple CPU cores + + Parameters + ---------- + b_array : numpy array [n_brightness_levels] + required brighntess levels. + alpha : TYPE + DESCRIPTION. + beta : TYPE + DESCRIPTION. + + Returns + ------- + sigmas : numpy array [n_brightness_levels] + Simulated std + diffs : numpy array [n_brightness_levels] + simulated differences + + """ if multiprocessing.cpu_count() > 1: N_CPU = multiprocessing.cpu_count()-1 else: @@ -70,12 +117,35 @@ def regular_MC(b_array, alpha, beta): diffs[b_] = result[0] sigmas[b_] = result[1] - + pool.close() return sigmas, diffs def interp_MC(b_array, sigma_min, sigma_max, diff_min, diff_max): + """ + Interpolates the missing values for diff and and sigma, based on their + upper and lower bounds. + + Parameters + ---------- + b_array : TYPE + DESCRIPTION. + sigma_min : TYPE + DESCRIPTION. + sigma_max : TYPE + DESCRIPTION. + diff_min : TYPE + DESCRIPTION. + diff_max : TYPE + DESCRIPTION. + + Returns + ------- + TYPE + DESCRIPTION. + + """ norm_b = (b_array - b_array[0])/(b_array[-1] - b_array[0]) sigmas_sq_lin = norm_b * (sigma_max**2 - sigma_min**2) + sigma_min**2 @@ -84,12 +154,32 @@ def interp_MC(b_array, return np.sqrt(sigmas_sq_lin[1:-1]), np.sqrt(diffs_sq_lin[1:-1]) -def run_fast_MC(alpha, beta): - t0 = time.perf_counter() - print("Beginning...") - xmin, xmax = get_non_linearity_bound(alpha, beta, tol) +def run_fast_MC(alpha, beta): + """ + Computes diff and sigmas for 1001 values of brighntess, for fixed alpha and + beta. - # these 2 inddexes define the points for which the linear model is valid. + To go faster, bound for which the probability of clipping is + neglectible are estimated (xmin, xmax). A regular MC is applied to brightness + outside of this range, and the rest of the points are interpolates linearly. + + Parameters + ---------- + alpha : TYPE + DESCRIPTION. + beta : TYPE + DESCRIPTION. + + Returns + ------- + None. + + """ + # t0 = time.perf_counter() + print("Estimating noise curves ...") + xmin, xmax = get_non_linearity_bound(alpha, beta, TOL) + + # these 2 indexes define the points for which the linear model is valid. # The regular Monte carlo is applied to these points, that are then used as reference # for the linear model imin = int(np.ceil(xmin * n_brightness_levels)) + 1 @@ -114,26 +204,27 @@ def run_fast_MC(alpha, beta): diffs[imin], diffs[imax]) sigmas[imin:imax+1] = sigmas_l diffs[imin:imax+1] = diffs_l - print('Finished! Estimated 1 noise curve in {:2f} sec'.format(time.perf_counter() - t0)) + # print('Finished! Estimated 1 noise curve in {:2f} sec'.format(time.perf_counter() - t0)) - plt.figure("sigma") - plt.plot(np.linspace(0, 1, n_brightness_levels+1), sigmas) - plt.grid() - plt.axline((xmin, 0), (xmin, 1e-5), c='k') - plt.axline((xmax, 0), (xmax, 1e-5), c='k') - plt.ylabel('Sigma') - plt.xlabel('brightness') + # plt.figure("sigma") + # plt.plot(np.linspace(0, 1, n_brightness_levels+1), sigmas) + # plt.grid() + # plt.axline((xmin, 0), (xmin, 1e-5), c='k') + # plt.axline((xmax, 0), (xmax, 1e-5), c='k') + # plt.ylabel('Sigma') + # plt.xlabel('brightness') - plt.figure("Diff") - plt.plot(np.linspace(0, 1, n_brightness_levels+1), diffs) - plt.grid() - plt.axline((xmin, 0), (xmin, 1e-5), c='k') - plt.axline((xmax, 0), (xmax, 1e-5), c='k') - plt.ylabel('Diff') - plt.xlabel('brightness') + # plt.figure("Diff") + # plt.plot(np.linspace(0, 1, n_brightness_levels+1), diffs) + # plt.grid() + # plt.axline((xmin, 0), (xmin, 1e-5), c='k') + # plt.axline((xmax, 0), (xmax, 1e-5), c='k') + # plt.ylabel('Diff') + # plt.xlabel('brightness') + return sigmas, diffs if __name__ == "__main__": diff --git a/handheld_super_resolution/super_resolution.py b/handheld_super_resolution/super_resolution.py index 4204c72..d74f311 100644 --- a/handheld_super_resolution/super_resolution.py +++ b/handheld_super_resolution/super_resolution.py @@ -15,30 +15,28 @@ """ import os -import glob import time from pathlib import Path import numpy as np from numba import cuda -import exifread import rawpy -from . import raw2rgb -from .utils import getTime, DEFAULT_NUMPY_FLOAT_TYPE, divide, add, round_iso from .utils_image import compute_grey_images, frame_count_denoising_gauss, frame_count_denoising_median -from .utils_dng import save_as_dng -from .merge import merge, merge_ref -from .kernels import estimate_kernels +from .utils import getTime, DEFAULT_NUMPY_FLOAT_TYPE, divide, add, round_iso from .block_matching import init_block_matching, align_image_block_matching -from .ICA import ICA_optical_flow, init_ICA -from .robustness import init_robustness, compute_robustness from .params import check_params_validity, get_params, merge_params +from .robustness import init_robustness, compute_robustness +from .utils_dng import save_as_dng, load_dng_burst +from .ICA import ICA_optical_flow, init_ICA +from .fast_monte_carlo import run_fast_MC +from .kernels import estimate_kernels +from .merge import merge, merge_ref +from . import raw2rgb NOISE_MODEL_PATH = Path(os.path.dirname(__file__)).parent / 'data' - def main(ref_img, comp_imgs, options, params): """ This is the implementation of Alg. 1: HandheldBurstSuperResolution. @@ -76,7 +74,7 @@ def main(ref_img, comp_imgs, options, params): accumulate_r = params['accumulated robustness denoiser']['on'] - #___ Moving to GPU + #### Moving to GPU cuda_ref_img = cuda.to_device(ref_img) cuda.synchronize() @@ -85,7 +83,7 @@ def main(ref_img, comp_imgs, options, params): t1 = time.perf_counter() - #___ Raw to grey + #### Raw to grey grey_method = params['grey method'] if bayer_mode : @@ -96,7 +94,7 @@ def main(ref_img, comp_imgs, options, params): else: cuda_ref_grey = cuda_ref_img - #___ Block Matching + #### Block Matching if verbose_2 : cuda.synchronize() current_time = time.perf_counter() @@ -109,7 +107,7 @@ def main(ref_img, comp_imgs, options, params): current_time = getTime(current_time, 'Block Matching initialised (Total)') - #___ ICA : compute grad and hessian + #### ICA : compute grad and hessian if verbose_2 : cuda.synchronize() current_time = time.perf_counter() @@ -122,7 +120,7 @@ def main(ref_img, comp_imgs, options, params): current_time = getTime(current_time, 'ICA initialised (Total)') - #___ Local stats estimation + #### Local stats estimation if verbose_2: cuda.synchronize() current_time = time.perf_counter() @@ -158,13 +156,13 @@ def main(ref_img, comp_imgs, options, params): print("\nProcessing image {} ---------\n".format(im_id+1)) im_time = time.perf_counter() - #___ Moving to GPU + #### Moving to GPU cuda_img = cuda.to_device(comp_imgs[im_id]) if verbose_3 : cuda.synchronize() current_time = getTime(im_time, 'Arrays moved to GPU') - #___ Compute Grey Images + #### Compute Grey Images if bayer_mode: cuda_im_grey = compute_grey_images(comp_imgs[im_id], grey_method) if verbose_3 : @@ -173,7 +171,7 @@ def main(ref_img, comp_imgs, options, params): else: cuda_im_grey = cuda_img - #___ Block Matching + #### Block Matching if verbose_2 : cuda.synchronize() current_time = time.perf_counter() @@ -186,7 +184,7 @@ def main(ref_img, comp_imgs, options, params): current_time = getTime(current_time, 'Block Matching (Total)') - #___ ICA + #### ICA if verbose_2 : cuda.synchronize() current_time = time.perf_counter() @@ -203,7 +201,7 @@ def main(ref_img, comp_imgs, options, params): current_time = getTime(current_time, 'Image aligned using ICA (Total)') - #___ Robustness + #### Robustness if verbose_2 : cuda.synchronize() current_time = time.perf_counter() @@ -219,7 +217,7 @@ def main(ref_img, comp_imgs, options, params): current_time = getTime(current_time, 'Robustness estimated (Total)') - #___ Kernel estimation + #### Kernel estimation if verbose_2 : cuda.synchronize() current_time = time.perf_counter() @@ -232,7 +230,7 @@ def main(ref_img, comp_imgs, options, params): current_time = getTime(current_time, 'Kernels estimated (Total)') - #___ Merging + #### Merging if verbose_2 : current_time = time.perf_counter() print('\nAccumulating Image') @@ -250,7 +248,7 @@ def main(ref_img, comp_imgs, options, params): if debug_mode : debug_dict['robustness'].append(cuda_robustness.copy_to_host()) - #___ Ref kernel estimation + #### Ref kernel estimation if verbose_2 : cuda.synchronize() current_time = time.perf_counter() @@ -262,7 +260,7 @@ def main(ref_img, comp_imgs, options, params): cuda.synchronize() current_time = getTime(current_time, 'Kernels estimated (Total)') - #___ Merge ref + #### Merge ref if verbose_2 : cuda.synchronize() print('\nAccumulating ref Img') @@ -323,103 +321,42 @@ def process(burst_path, options=None, custom_params=None): options['verbose'] >= 2) params = {} - ref_id = 0 #TODO Select ref id based on HDR+ method - - raw_comp = [] - - # This ensures that burst_path is a Path object - burst_path = Path(burst_path) + # reading image stack + ref_raw, raw_comp, ISO, tags, CFA, xyz2cam, ref_path = load_dng_burst(burst_path) - - # Get the list of raw images in the burst path - raw_path_list = glob.glob(os.path.join(burst_path.as_posix(), '*.dng')) - assert len(raw_path_list) != 0, 'At least one raw .dng file must be present in the burst folder.' - # Read the raw bayer data from the DNG files - for index, raw_path in enumerate(raw_path_list): - with rawpy.imread(raw_path) as rawObject: - if index != ref_id : - - raw_comp.append(rawObject.raw_image.copy()) # copy otherwise image data is lost when the rawpy object is closed - raw_comp = np.array(raw_comp) - - # Reference image selection and metadata - raw = rawpy.imread(raw_path_list[ref_id]) - ref_raw = raw.raw_image.copy() - xyz2cam = raw2rgb.get_xyz2cam_from_exif(raw_path_list[ref_id]) - - - # reading exifs for white level, black leve and CFA - with open(raw_path_list[ref_id], 'rb') as raw_file: - tags = exifread.process_file(raw_file) - - - white_level = tags['Image Tag 0xC61D'].values[0] # there is only one white level - - black_levels = tags['Image BlackLevel'] - if isinstance(black_levels.values[0], int): - black_levels = np.array(black_levels.values) - else: # Sometimes this tag is a fraction object for some reason. It seems that black levels are all integers anyway - black_levels = np.array([int(x.decimal()) for x in black_levels.values]) - - white_balance = raw.camera_whitebalance - - CFA = tags['Image CFAPattern'] - CFA = np.array(list(CFA.values)).reshape(2,2) - - if 'EXIF ISOSpeedRatings' in tags.keys(): - ISO = int(str(tags['EXIF ISOSpeedRatings'])) - elif 'Image ISOSpeedRatings' in tags.keys(): - ISO = int(str(tags['Image ISOSpeedRatings'])) + # if the algorithm had to be run on a specific sensor, + # the precise values of alpha and beta could be used instead + if 'mode' in custom_params and custom_params['mode'] == 'grey': + alpha = tags['Image Tag 0xC761'].values[0][0] + beta = tags['Image Tag 0xC761'].values[1][0] else: - raise AttributeError('ISO value could not be found in both EXIF and Image type.') - - # Clipping ISO to 100 from below - ISO = max(100, ISO) - ISO = min(3200, ISO) + # Averaging RGB noise values + ## IMPORTANT NOTE : the noise model exif already are NOT for nominal ISO 100 + ## But are already scaled for the image ISO. + alpha = sum([x[0] for x in tags['Image Tag 0xC761'].values[::2]])/3 + beta = sum([x[0] for x in tags['Image Tag 0xC761'].values[1::2]])/3 + + #### Packing noise model related to picture ISO + # curve_iso = round_iso(ISO) # Rounds non standart ISO to regular ISO (100, 200, 400, ...) + # std_noise_model_label = 'noise_model_std_ISO_{}'.format(curve_iso) + # diff_noise_model_label = 'noise_model_diff_ISO_{}'.format(curve_iso) + # std_noise_model_path = (NOISE_MODEL_PATH / std_noise_model_label).with_suffix('.npy') + # diff_noise_model_path = (NOISE_MODEL_PATH / diff_noise_model_label).with_suffix('.npy') - # Packing noise model related to picture ISO - curve_iso = round_iso(ISO) # Rounds non standart ISO to regular ISO (100, 200, 400, ...) - std_noise_model_label = 'noise_model_std_ISO_{}'.format(curve_iso) - diff_noise_model_label = 'noise_model_diff_ISO_{}'.format(curve_iso) - std_noise_model_path = (NOISE_MODEL_PATH / std_noise_model_label).with_suffix('.npy') - diff_noise_model_path = (NOISE_MODEL_PATH / diff_noise_model_label).with_suffix('.npy') + # std_curve = np.load(std_noise_model_path) + # diff_curve = np.load(diff_noise_model_path) - std_curve = np.load(std_noise_model_path) - diff_curve = np.load(diff_noise_model_path) + ## Use this to compute noise curves on the fly + std_curve, diff_curve = run_fast_MC(alpha, beta) if verbose_2: currentTime = getTime(currentTime, ' -- Read raw files') - - if np.issubdtype(type(ref_raw[0,0]), np.integer): - ## Here do black and white level correction and white balance processing for all image in comp_images - ## Each image in comp_images should be between 0 and 1. - ## ref_raw is a (H,W) array - ref_raw = ref_raw.astype(DEFAULT_NUMPY_FLOAT_TYPE) - for i in range(2): - for j in range(2): - channel = CFA[i, j] - ref_raw[i::2, j::2] = (ref_raw[i::2, j::2] - black_levels[channel]) / (white_level - black_levels[channel]) - ref_raw[i::2, j::2] *= white_balance[channel] / white_balance[1] - - - ref_raw = np.clip(ref_raw, 0.0, 1.0) - ## The division by the green WB value is important because WB may come with integer coefficients instead - - if np.issubdtype(type(raw_comp[0,0,0]), np.integer): - raw_comp = raw_comp.astype(DEFAULT_NUMPY_FLOAT_TYPE) - ## raw_comp is a (N, H,W) array - for i in range(2): - for j in range(2): - channel = channel = CFA[i, j] - raw_comp[:, i::2, j::2] = (raw_comp[:, i::2, j::2] - black_levels[channel]) / (white_level - black_levels[channel]) - raw_comp[:, i::2, j::2] *= white_balance[channel] / white_balance[1] - raw_comp = np.clip(raw_comp, 0., 1.) - - #___ Estimating ref image SNR + + #### Estimating ref image SNR brightness = np.mean(ref_raw) id_noise = round(1000*brightness) @@ -435,7 +372,7 @@ def process(burst_path, options=None, custom_params=None): SNR_params = get_params(SNR) - #__ Merging params dictionnaries + #### Merging params dictionnaries # checking (just in case !) check_params_validity(SNR_params, ref_raw.shape) @@ -444,21 +381,10 @@ def process(burst_path, options=None, custom_params=None): params = merge_params(dominant=custom_params, recessive=SNR_params) check_params_validity(params, ref_raw.shape) - #__ adding metadatas to dict + #### adding metadatas to dict if not 'noise' in params['merging'].keys(): params['merging']['noise'] = {} - # if the algorithm had to be run on a specific sensor, - # the precise values of alpha and beta could be used instead - if params['mode'] == 'grey': - alpha = tags['Image Tag 0xC761'].values[0][0] - beta = tags['Image Tag 0xC761'].values[1][0] - else: - # Averaging RGB noise values - ## IMPORTAN0T NOTE : the noise model exif already are NOT for nominal ISO 100 - ## But are already scaled for the image ISO. - alpha = sum([x[0] for x in tags['Image Tag 0xC761'].values[::2]])/3 - beta = sum([x[0] for x in tags['Image Tag 0xC761'].values[1::2]])/3 params['merging']['noise']['alpha'] = alpha params['merging']['noise']['beta'] = beta @@ -513,24 +439,16 @@ def process(burst_path, options=None, custom_params=None): params['merging']['accumulated robustness denoiser'] = params['accumulated robustness denoiser']['merge'] else: params['merging']['accumulated robustness denoiser'] = {'on' : False} - - + - #___ Running the handheld pipeline + #### Running the handheld pipeline handheld_output, debug_dict = main(ref_raw.astype(DEFAULT_NUMPY_FLOAT_TYPE), raw_comp.astype(DEFAULT_NUMPY_FLOAT_TYPE), options, params) - #___ Saving output as dng if required - if 'dng_outpath' in options.keys(): - outpath = Path(options["dng_outpath"]) - if verbose_1: - print('Saving output to {}'.format(outpath.with_suffix('.dng').as_posix())) - save_as_dng(handheld_output.copy_to_host(), raw_path_list[ref_id], outpath) - - #___ Performing frame count aware denoising if enabled + #### Performing frame count aware denoising if enabled median_params = params['accumulated robustness denoiser']['median'] gauss_params = params['accumulated robustness denoiser']['gauss'] @@ -557,12 +475,13 @@ def process(burst_path, options=None, custom_params=None): gauss_params) - #___ post processing + #### post processing if post_processing_enabled: if verbose_2: print('-- Post processing image') - + + raw = rawpy.imread(ref_path) output_image = raw2rgb.postprocess(raw, handheld_output.copy_to_host(), params_pp['do color correction'], params_pp['do tonemapping'], @@ -575,7 +494,7 @@ def process(burst_path, options=None, custom_params=None): else: output_image = handheld_output.copy_to_host() - #__ return + if params['debug']: return output_image, debug_dict diff --git a/handheld_super_resolution/utils_dng.py b/handheld_super_resolution/utils_dng.py index a076b6c..2930a13 100644 --- a/handheld_super_resolution/utils_dng.py +++ b/handheld_super_resolution/utils_dng.py @@ -6,14 +6,137 @@ """ import os +import glob import numpy as np +from pathlib import Path import exifread import rawpy import imageio -EXIFTOOL_PATH = 'C:/Users/jamyl/Downloads/exiftool-12.44/exiftool.exe' -DNG_VALIDATE_PATH = 'C:/Users/jamyl/Downloads/dng_sdk_1_6/dng_sdk_1_6/dng_sdk/targets/win/debug64_x64/dng_validate.exe' +from . import raw2rgb +from .utils import DEFAULT_NUMPY_FLOAT_TYPE + +# Paths of exiftool and dng validate. Only necessary to output dng. +EXIFTOOL_PATH = 'path/to/exiftool.exe' +DNG_VALIDATE_PATH = 'path/to/dng_validate.exe' + + +def load_dng_burst(burst_path): + """ + Loads a dng burst into numpy arrays, and their exif tags. + + Parameters + ---------- + burst_path : Path or str + Path of the folder containing the .dngs + + Returns + ------- + ref_raw : numpy Array[H, W] + Reference frame + raw_comp : numpy Array[n, H, W] + Stack of non-reference frame + ISO : int + Clipped ISO (between 100 and 3600) + tags : dict + Tags of the reference frame + CFA : numpy array [2, 2] + Bayer pattern of the stack + xyz2cam : TYPE + DESCRIPTION. + TYPE + DESCRIPTION. + + """ + ref_id = 0 # TODO Select ref id based on HDR+ method + raw_comp = [] + + # This ensures that burst_path is a Path object + burst_path = Path(burst_path) + + + #### Read dng as numpy arrays + # Get the list of raw images in the burst path + raw_path_list = glob.glob(os.path.join(burst_path.as_posix(), '*.dng')) + assert len(raw_path_list) != 0, 'At least one raw .dng file must be present in the burst folder.' + # Read the raw bayer data from the DNG files + for index, raw_path in enumerate(raw_path_list): + with rawpy.imread(raw_path) as rawObject: + if index != ref_id: + + raw_comp.append(rawObject.raw_image.copy()) # copy otherwise image data is lost when the rawpy object is closed + raw_comp = np.array(raw_comp) + + # Reference image selection and metadata + raw = rawpy.imread(raw_path_list[ref_id]) + ref_raw = raw.raw_image.copy() + + + + + #### Reading tags of the reference image + xyz2cam = raw2rgb.get_xyz2cam_from_exif(raw_path_list[ref_id]) + + # reading exifs for white level, black leve and CFA + with open(raw_path_list[ref_id], 'rb') as raw_file: + tags = exifread.process_file(raw_file) + + white_level = tags['Image Tag 0xC61D'].values[0] # there is only one white level + + black_levels = tags['Image BlackLevel'] + if isinstance(black_levels.values[0], int): + black_levels = np.array(black_levels.values) + else: # Sometimes this tag is a fraction object for some reason. It seems that black levels are all integers anyway + black_levels = np.array([int(x.decimal()) for x in black_levels.values]) + + white_balance = raw.camera_whitebalance + + CFA = tags['Image CFAPattern'] + CFA = np.array(list(CFA.values)).reshape(2, 2) + + if 'EXIF ISOSpeedRatings' in tags.keys(): + ISO = int(str(tags['EXIF ISOSpeedRatings'])) + elif 'Image ISOSpeedRatings' in tags.keys(): + ISO = int(str(tags['Image ISOSpeedRatings'])) + else: + raise AttributeError('ISO value could not be found in both EXIF and Image type.') + + # Clipping ISO to 100 from below + ISO = max(100, ISO) + ISO = min(3200, ISO) + + + + + #### Performing whitebalance and normalizing into 0, 1 + + if np.issubdtype(type(ref_raw[0, 0]), np.integer): + # Here do black and white level correction and white balance processing for all image in comp_images + # Each image in comp_images should be between 0 and 1. + # ref_raw is a (H,W) array + ref_raw = ref_raw.astype(DEFAULT_NUMPY_FLOAT_TYPE) + for i in range(2): + for j in range(2): + channel = CFA[i, j] + ref_raw[i::2, j::2] = (ref_raw[i::2, j::2] - black_levels[channel]) / (white_level - black_levels[channel]) + ref_raw[i::2, j::2] *= white_balance[channel] / white_balance[1] + + ref_raw = np.clip(ref_raw, 0.0, 1.0) + # The division by the green WB value is important because WB may come with integer coefficients instead + + if np.issubdtype(type(raw_comp[0, 0, 0]), np.integer): + raw_comp = raw_comp.astype(DEFAULT_NUMPY_FLOAT_TYPE) + # raw_comp is a (N, H,W) array + for i in range(2): + for j in range(2): + channel = channel = CFA[i, j] + raw_comp[:, i::2, j::2] = (raw_comp[:, i::2, j::2] - black_levels[channel]) / (white_level - black_levels[channel]) + raw_comp[:, i::2, j::2] *= white_balance[channel] / white_balance[1] + raw_comp = np.clip(raw_comp, 0., 1.) + + return ref_raw, raw_comp, ISO, tags, CFA, xyz2cam, raw_path_list[ref_id] + def save_as_dng(np_img, ref_dng_path, outpath): ''' @@ -21,15 +144,15 @@ def save_as_dng(np_img, ref_dng_path, outpath): The image is first saved as 16bits tiff, then the extension is swapped to .dng. The metadata are then overwritten using a reference dng, and the final dng is built using dng_validate. - + Requires : - dng_validate (can be found in dng sdk): https://helpx.adobe.com/camera-raw/digital-negative.html#dng_sdk_download - + - exiftool https://exiftool.org/ - - + + Based on : https://github.com/gluijk/dng-from-tiff/blob/main/dngmaker.bat https://github.com/antonwolf/dng_stacker/blob/master/dng_stacker.bat @@ -49,45 +172,41 @@ def save_as_dng(np_img, ref_dng_path, outpath): ''' assert np_img.shape[-1] == 3 - - - np_int_img = np.copy(np_img) # copying to avoid inplace-overwritting - # Undo White balance and black level - ## get tags + + np_int_img = np.copy(np_img) # copying to avoid inplace-overwritting + #### Undo White balance and black level + # get tags with open(ref_dng_path, 'rb') as raw_file: tags = exifread.process_file(raw_file) - white_level = tags['Image Tag 0xC61D'].values[0] # there is only one white level - - black_levels = tags['Image BlackLevel'] + black_levels = tags['Image BlackLevel'] if isinstance(black_levels.values[0], int): black_levels = np.array(black_levels.values) - else: # Sometimes this tag is a fraction object for some reason. It seems that black levels are all integers anyway + else: # Sometimes this tag is a fraction object for some reason. It seems that black levels are all integers anyway black_levels = np.array([int(x.decimal()) for x in black_levels.values]) - + raw = rawpy.imread(ref_dng_path) white_balance = raw.camera_whitebalance - - ## Reverse WB - + + # Reverse WB + new_white_level = 2**16 - 1 + for c in range(3): np_int_img[:, :, c] /= white_balance[c] / white_balance[1] - np_int_img[:, :, c] = np_int_img[:, :, c] * (white_level - black_levels[c]) + black_levels[c] + np_int_img[:, :, c] = np_int_img[:, :, c] * (new_white_level - black_levels[c]) + black_levels[c] + + np_int_img = np.clip(np_int_img, 0, 2**16 - 1).astype(np.uint16) - np_int_img = np.clip(np_int_img, 0, 2**16-1).astype(np.uint16) - - # Saving the image as 16 bits RGB tiff + #### Saving the image as 16 bits RGB tiff save_as_tiff(np_int_img, outpath) - - tmp_path = outpath.parent/'tmp.dng' - + + tmp_path = outpath.parent / 'tmp.dng' + # Deleting tmp.dng if it is already existing if os.path.exists(tmp_path): os.remove(tmp_path) - - - - # Overwritting the tiff tags with dng tags, and replacing the .tif extension + + #### Overwritting the tiff tags with dng tags, and replacing the .tif extension # by .dng cmd = ''' {} -n\ @@ -114,10 +233,10 @@ def save_as_dng(np_img, ref_dng_path, outpath): "-IFD0:OpcodeList3 2: - print(' WARNING: Since the optics and the integration on the sensor limit the aliasing, do not expect more details than that obtained at x2 (refer to our paper and the original publication).') -print('') -if args.R_on: - print(' Robustness: enabled') + print(' Super-resolution mode.') + if args.scale > 2: + print(' WARNING: Since the optics and the integration on the sensor limit the aliasing, do not expect more details than that obtained at x2 (refer to our paper and the original publication).') + print('') + if args.R_on: + print(' Robustness: enabled') + print(' -------------------------') + print(' t: %1.2f' % args.t) + print(' s1: %1.2f' % args.s1) + print(' s2: %1.2f' % args.s2) + print(' Mt: %1.2f' % args.Mt) + if args.R_denoising_on: + print(' Robustness denoising: enabled') + else: + print(' Robustness denoising: disabled') + print(' ' ) + else: + print(' Robustness: disabled') + print(' ' ) + + print(' Alignment:') + print(' -------------------------') + print(' ICA Iterations: %d' % args.ICA_iter) + print('') + print(' Fusion:') print(' -------------------------') - print(' t: %1.2f' % args.t) - print(' s1: %1.2f' % args.s1) - print(' s2: %1.2f' % args.s2) - print(' Mt: %1.2f' % args.Mt) - if args.R_denoising_on: - print(' Robustness denoising: enabled') + print(' Kernel shape: %s' % args.kernel_shape) + print(' k_stretch: %1.2f' % args.k_stretch) + print(' k_shrink: %1.2f' % args.k_shrink) + if args.k_detail is not None: + print(' k_detail: %1.2f' % args.k_detail) else: - print(' Robustness denoising: disabled') - print(' ' ) -else: - print(' Robustness: disabled') - print(' ' ) - -print(' Alignment:') -print(' -------------------------') -print(' ICA Iterations: %d' % args.ICA_iter) -print('') -print(' Fusion:') -print(' -------------------------') -print(' Kernel shape: %s' % args.kernel_shape) -print(' k_stretch: %1.2f' % args.k_stretch) -print(' k_shrink: %1.2f' % args.k_shrink) -if args.k_detail is not None: - print(' k_detail: %1.2f' % args.k_detail) -else: - print(' k_detail: SNR based' ) -if args.k_denoise is not None: - print(' k_denoise: %1.2f' % args.k_denoise) -else: - print(' k_denoise: SNR based' ) -print('') - - - - - -#### Handheld #### -print('Processing with handheld super-resolution') -options = {'verbose' : args.verbose} -params={"scale": args.scale, - "merging":{"kernel":args.kernel_shape}, - "robustness":{"on":args.R_on}, - "kanade": {"tuning": {"kanadeIter":args.ICA_iter}} - } - - -params['robustness']['tuning'] = {'t' : args.t, - 's1' : args.s1, - 's2' : args.s2, - 'Mt' : args.Mt, - } - -if args.k_detail is not None: - params['merging']['tuning']['k_detail'] = args.k_detail -if args.k_denoise is not None: - params['merging']['tuning']['k_denoise'] = args.k_denoise - -params['merging'] = {'tuning' : {'k_stretch' : args.k_stretch, - 'k_shrink' : args.k_shrink - }} -params['accumulated robustness denoiser'] = {'on': args.R_denoising_on} - -params['post processing'] = {'on':args.post_process, - 'do sharpening' : args.do_sharpening, - 'do tonemapping':args.do_tonemapping, - 'do gamma' : args.do_gamma, - 'do devignette' : False, - 'do color correction': args.do_color_correction, + print(' k_detail: SNR based' ) + if args.k_denoise is not None: + print(' k_denoise: %1.2f' % args.k_denoise) + else: + print(' k_denoise: SNR based' ) + print('') + + + + + + #### Handheld #### + print('Processing with handheld super-resolution') + options = {'verbose' : args.verbose} + params={"scale": args.scale, + "merging":{"kernel":args.kernel_shape}, + "robustness":{"on":args.R_on}, + "kanade": {"tuning": {"kanadeIter":args.ICA_iter}} + } + + + params['robustness']['tuning'] = {'t' : args.t, + 's1' : args.s1, + 's2' : args.s2, + 'Mt' : args.Mt, + } + + if args.k_detail is not None: + params['merging']['tuning']['k_detail'] = args.k_detail + if args.k_denoise is not None: + params['merging']['tuning']['k_denoise'] = args.k_denoise + + params['merging'] = {'tuning' : {'k_stretch' : args.k_stretch, + 'k_shrink' : args.k_shrink + }} + params['accumulated robustness denoiser'] = {'on': args.R_denoising_on} + + outpath = Path(args.outpath) + # disabling post processing for dng outputs + if outpath.suffix == '.dng': + args.post_process = False + + + params['post processing'] = {'on':args.post_process, + 'do sharpening' : args.do_sharpening, + 'do tonemapping':args.do_tonemapping, + 'do gamma' : args.do_gamma, + 'do devignette' : False, + 'do color correction': args.do_color_correction, + + 'sharpening' : {'radius': args.radius, + 'ammount': args.amount} + } + + + handheld_output = process(args.impath, options, params) + handheld_output = np.nan_to_num(handheld_output) + handheld_output = np.clip(handheld_output, 0, 1) + + + # define a faster imsave for large png images + def imsave(fname, rgb_8bit_data): + return cv2.imwrite(fname, cv2.cvtColor(rgb_8bit_data, cv2.COLOR_RGB2BGR )) + + + #### Save images #### + + if outpath.suffix == '.dng': + if options['verbose'] >=1 : + print('Saving output to {}'.format(outpath.with_suffix('.dng').as_posix())) + ref_img_path = glob.glob(os.path.join(args.impath, '*.dng'))[0] + save_as_dng(handheld_output, ref_img_path, outpath) - 'sharpening' : {'radius': args.radius, - 'ammount': args.amount} - } - - -handheld_output = process(args.impath, options, params) -handheld_output = np.nan_to_num(handheld_output) -handheld_output = np.clip(handheld_output, 0, 1) - - -# define a faster imsave for large png images -def imsave(fname, rgb_8bit_data): - return cv2.imwrite(fname, cv2.cvtColor(rgb_8bit_data, cv2.COLOR_RGB2BGR )) - - -#### Save images #### -imsave(args.outpath, img_as_ubyte(handheld_output)) - -# Robustness map -print('done') \ No newline at end of file + else: + imsave(args.outpath, img_as_ubyte(handheld_output)) \ No newline at end of file From 32ede7e24dac1517eec066f1ba1d30602ee3af54 Mon Sep 17 00:00:00 2001 From: Jamy Date: Mon, 27 Mar 2023 15:00:00 +0200 Subject: [PATCH 4/4] Update README.md --- README.md | 34 ++++++++++++++++++++++++++++++++++ 1 file changed, 34 insertions(+) diff --git a/README.md b/README.md index 72fd2b5..639c7f2 100644 --- a/README.md +++ b/README.md @@ -65,6 +65,40 @@ The core of the algorithm is in `handheld_super_resolution.process(burst_path, o To obtain the bursts used in the publication, please download the latest release of the repo. It contains the code and two raw bursts of respectively 13 images from [[Bhat et al., ICCV21]](https://arxiv.org/abs/2108.08286) and 20 images from [[Lecouat et al., SIGGRAPH22]](https://arxiv.org/abs/2207.14671). Otherwise specify the path to any burst of raw images, e.g., `*.dng`, `*.ARW` or `*.CR2` for instance. The result is found in the `./results/` folder. Remember that if you have activated the post-processing flag, the predicted image will be further tone-mapped and sharpened. Deactivate it if you want to plug in your own ISP. +## Generating DNG + +Saving images as DNG requires extra-dependencies and steps: +1. Install imageio using pip: + + ```bash + pip install imageio + ``` + +2. Install exiftool. You can download the appropriate version for your operating system from the [exiftool website](https://exiftool.org/). + +3. Download the DNG SDK from the [Adobe website](https://helpx.adobe.com/camera-raw/digital-negative.html#dng_sdk_download). + +4. Follow the instructions in the `DNG_ReadMe.txt` file included in the downloaded DNG SDK to complete the installation. + +5. Set the paths to the exiftool and DNG validate executables in the `utils_dng.py` script. Open the `utils_dng.py` file located in the project directory, and modify the values of the `EXIFTOOL_PATH` and `DNG_VALIDATE_PATH` variables to match the paths to the exiftool and DNG validate executables on your system. + + ```python + # utils_dng.py + + # Set the path to the exiftool executable + EXIFTOOL_PATH = "/path/to/exiftool" + + # Set the path to the DNG validate executable + DNG_VALIDATE_PATH = "/path/to/dng_validate" + ``` + + Note that you may need to adjust the paths based on the location where you installed exiftool and DNG validate on your system. + +That's it! Once you've completed these steps you should be able to save the output as DNG, for example by running: +``` +python run_handheld.py --impath test_burst --outpath output.dng +``` + ## Citation If this code or the implementation details of the companion IPOL publication are of any help, please cite our work: ```BibTex