100 lines
4.5 KiB
Python
100 lines
4.5 KiB
Python
# Notes: https://gitea.lemnoslife.com/Benjamin_Loison/Robust_image_source_identification_on_modern_smartphones/issues/25
|
|
|
|
import os
|
|
from PIL import Image
|
|
import numpy as np
|
|
import matplotlib.pyplot as plt
|
|
import sys
|
|
|
|
sys.path.insert(0, '../../algorithms/image_utils/')
|
|
|
|
from image_utils import showImageWithMatplotlib, randomGaussianImage, toPilImage, getPrnuShownAsSuch
|
|
|
|
sys.path.insert(0, '../../algorithms/context_adaptive_interpolator/')
|
|
|
|
from context_adaptive_interpolator import contextAdaptiveInterpolator
|
|
|
|
sys.path.insert(0, '../../algorithms/distance/')
|
|
|
|
from rms_diff import rmsDiffNumpy
|
|
|
|
from skimage.restoration import denoise_tv_chambolle
|
|
|
|
datasetPath = 'no_noise_images'
|
|
# Note that contrarily to `datasets/fake/`, here we do not have images being Gaussian with `scale` `1` but actual images with pixel values between 0 and 255.
|
|
# In addition to the range difference, note that the distribution in the first set of images was a Gaussian and here is very different and specific.
|
|
PRNU_FACTOR = 0.01
|
|
NOISE_FACTOR = 0.1
|
|
|
|
np.random.seed(0)
|
|
|
|
SPLIT_N_X_N_S = [1, 2, 4]
|
|
|
|
# len(SPLIT_N_X_N_S)
|
|
fig, axes = plt.subplots(2, 4)
|
|
fig.suptitle('PRNU estimation with different number of images having Gaussian noise and Gaussian noised PRNU')
|
|
|
|
for splitNXNIndex, splitNXN in enumerate(SPLIT_N_X_N_S):
|
|
IMAGE_SIZE_SHAPE = [dimension // splitNXN for dimension in (704, 469)]
|
|
|
|
#prnuNpArray = 255 * randomGaussianImage(scale = PRNU_FACTOR, size = IMAGE_SIZE_SHAPE)
|
|
prnuNpArray = getPrnuShownAsSuch(IMAGE_SIZE_SHAPE, 255) * PRNU_FACTOR
|
|
|
|
def isIn256Range(x):
|
|
return 0 <= x and x <= 255
|
|
|
|
imagesPrnuEstimateNpArray = []
|
|
isFirstImage = True
|
|
|
|
for imageName in os.listdir(datasetPath):
|
|
if imageName.endswith('.png'):
|
|
imagePath = f'{datasetPath}/{imageName}'
|
|
imageWithoutPrnuPil = Image.open(imagePath).convert('F')
|
|
imageWithoutPrnuNpArray = np.array(imageWithoutPrnuPil)
|
|
|
|
m = IMAGE_SIZE_SHAPE[1]
|
|
n = IMAGE_SIZE_SHAPE[0]
|
|
|
|
imageWithoutPrnuNpArrayTiles = [imageWithoutPrnuNpArray[x : x + m, y : y + n] for x in range(0, imageWithoutPrnuNpArray.shape[0], m) for y in range(0, imageWithoutPrnuNpArray.shape[1], n)]
|
|
for imageWithoutPrnuNpArrayTile in imageWithoutPrnuNpArrayTiles:
|
|
#print(imageWithoutPrnuNpArrayTile.shape, tuple(IMAGE_SIZE_SHAPE[::-1]))
|
|
#if imageWithoutPrnuNpArrayTile.shape != tuple(IMAGE_SIZE_SHAPE[::-1]):
|
|
# continue
|
|
imageNoise = randomGaussianImage(scale = 255 * NOISE_FACTOR, size = imageWithoutPrnuNpArrayTile.shape)
|
|
imageWithPrnuNpArray = imageWithoutPrnuNpArrayTile + prnuNpArray + imageNoise
|
|
|
|
if splitNXNIndex == 0 and isFirstImage:
|
|
axis = axes[0]
|
|
|
|
axis[0].set_title('First image without noise')
|
|
axis[0].imshow(imageWithoutPrnuNpArrayTile)
|
|
|
|
axis[1].set_title('Actual Gaussian noised PRNU')
|
|
axis[1].imshow(prnuNpArray)
|
|
|
|
axis[2].set_title('F. i. with G. n.')
|
|
axis[2].imshow(imageWithoutPrnuNpArray + imageNoise)
|
|
|
|
axis[3].set_title('F. i. with G. n. and PRNU')
|
|
axis[3].imshow(imageWithoutPrnuNpArray + prnuNpArray + imageNoise)
|
|
isFirstImage = False
|
|
|
|
#assert all([isIn256Range(extreme) for extreme in [imageWithPrnuNpArray.max(), imageWithPrnuNpArray.min()]]), 'Adding the PRNU resulted in out of 256 bounds image'
|
|
imageWithPrnuPil = toPilImage(imageWithPrnuNpArray)
|
|
#imagePrnuEstimatePil = contextAdaptiveInterpolator(imageWithPrnuPil.load(), imageWithPrnuPil)
|
|
#imagePrnuEstimateNpArray = np.array(imagePrnuEstimatePil)
|
|
imagePrnuEstimateNpArray = imageWithPrnuNpArray - denoise_tv_chambolle(imageWithPrnuNpArray, weight=0.2, channel_axis=-1)
|
|
|
|
imagesPrnuEstimateNpArray += [imagePrnuEstimateNpArray]
|
|
|
|
cameraPrnuEstimateNpArray = np.array(imagesPrnuEstimateNpArray).mean(axis = 0)
|
|
rms = rmsDiffNumpy(cameraPrnuEstimateNpArray, prnuNpArray, True)
|
|
title = f'RMS with actual PRNU: {rmsDiffNumpy(cameraPrnuEstimateNpArray, prnuNpArray):.4f}\n(normalized RMS: {rmsDiffNumpy(cameraPrnuEstimateNpArray, prnuNpArray, True):.4f})'
|
|
axis = axes[1]
|
|
axis[splitNXNIndex].set_title(f'Number of images: {len(imagesPrnuEstimateNpArray)}\n{title}')
|
|
axis[splitNXNIndex].imshow(cameraPrnuEstimateNpArray)
|
|
|
|
axes[1][3].axis('off')
|
|
|
|
plt.tight_layout()
|
|
plt.show() |