From 7db83a10c9d12afd043e77df254d2b85654d4c8f Mon Sep 17 00:00:00 2001 From: Benjamin Loison Date: Thu, 23 May 2024 12:17:54 +0200 Subject: [PATCH] Update uncommited --- datasets/raise/fft/remove_period_patterns.py | 30 ++++++++++---------- 1 file changed, 15 insertions(+), 15 deletions(-) diff --git a/datasets/raise/fft/remove_period_patterns.py b/datasets/raise/fft/remove_period_patterns.py index d27863e..8f591e4 100644 --- a/datasets/raise/fft/remove_period_patterns.py +++ b/datasets/raise/fft/remove_period_patterns.py @@ -24,43 +24,43 @@ image = multipleColorsImage fft1 = fftpack.fftshift(fftpack.fft2(image)) -originalFft1 = fft1.copy() - def removePeriodicPatterns(fft1Part): # This example is intended to demonstrate how astropy.convolve and # scipy.convolve handle missing data, so we start by setting the brightest # pixels to NaN to simulate a "saturated" data set + # Unclear why the Gaussian kernel usage does not work as expected, even when just consider the *real* part. # See [Benjamin_Loison/astropy/issues/2](https://codeberg.org/Benjamin_Loison/astropy/issues/2). + fixedImagePart = fft1Part.copy() height, width = fft1Part.shape for x in range(width): - fft1Part[height // 2, x] = np.nan + #fft1Part[height // 2, x] = np.nan + fixedImagePart[height // 2, x] = np.mean([fft1Part[height // 2 - 1, x], fft1Part[height // 2 + 1, x]]) middleX = width // 2 RANGE = 1 for verticalLineX in [width // 4, width // 2, round(3 * width / 4)]: for y in range(height): for x in range(verticalLineX - RANGE, verticalLineX + RANGE + 1): - fft1Part[y, x] = np.nan + #fft1Part[y, x] = np.nan + fixedImagePart[y, x] = np.mean([fft1Part[y, x - RANGE - 1], fft1Part[y, x + RANGE + 1]]) # We smooth with a Gaussian kernel - kernel = Gaussian2DKernel(x_stddev = X_STDDEV) + #kernel = Gaussian2DKernel(x_stddev = X_STDDEV) # create a "fixed" image with NaNs replaced by interpolated values # See [Benjamin_Loison/astropy/issues/1](https://codeberg.org/Benjamin_Loison/astropy/issues/1). - fixedImagePart = interpolate_replace_nans(fft1Part, kernel) + #fixedImagePart = interpolate_replace_nans(fft1Part, kernel) return fixedImagePart -# TODO: are `.copy()` really necessary? -realFixedImage = removePeriodicPatterns(np.real(fft1).copy()) -imaginaryFixedImage = removePeriodicPatterns(np.imag(fft1).copy()) +realFixedImage = removePeriodicPatterns(np.real(fft1)) +imaginaryFixedImage = removePeriodicPatterns(np.imag(fft1)) fixedImage = realFixedImage + 1j * imaginaryFixedImage -#plt.imsave('second_image.png', np.log10(1 + abs(fixedImage))) figure, axes = plt.subplots(1, 2, sharex = True, sharey = True) plt.suptitle('Attenuating FFT significant lines') axes[0].set_title('Original FFT') -firstImage = np.log10(1 + abs(originalFft1)) +firstImage = np.log10(1 + abs(fft1)) secondImage = np.log10(1 + abs(fixedImage)) images = [firstImage, secondImage] vMin = np.min(images) @@ -77,18 +77,18 @@ def inverseFft(fft): ifft2 = np.real(fftpack.ifft2(fftpack.ifftshift(fft))) return ifft2 -invOriginalFft1 = inverseFft(originalFft1) +invFft1 = inverseFft(fft1) invFixedImage = inverseFft(fixedImage) -images = [invOriginalFft1, invFixedImage] +images = [invFft1, invFixedImage] vMin = np.min(images) vMax = np.max(images) plt.suptitle('Rafael 23/04/24 PRNU mean denoiser with periodic patterns attenuated') axes[0].set_title('Original PRNU') -axes[0].imshow(invOriginalFft1, vmin = vMin, vmax = vMax) +axes[0].imshow(invFft1, vmin = vMin, vmax = vMax) axes[1].set_title('PRNU with periodic patterns attenuated') axes[1].imshow(invFixedImage, vmin = vMin, vmax = vMax) axes[2].set_title('Difference between both left images') -differenceBetweenBothImages = invOriginalFft1 - invFixedImage +differenceBetweenBothImages = invFft1 - invFixedImage # `np.log10(1 + abs(differenceBetweenBothImages))` does not seem more interesting. axes[2].imshow(differenceBetweenBothImages) plt.tight_layout()