实验报告
实验题目：中值滤波取证与重采样取证
课程题目：多媒体安全技术
学号：5120180885
1、实验目的和要求
本次实验分两个部分：中值滤波实验和重采样实验
中值滤波：
任务1：根据中值滤波论文实现 公式(3), 和 ho_hat = ho / y ，并用于区分实验3.zip/img/中50张原图图像。（自行进行中值滤波处理）。画出类论文Figure.6的图。
任务2：根据论文实现公式(4), q = h0 / h1 , 同理区分50张原图。画出Figure.6
任务3：根据论文实现公式(5),(6)。 q_hat = median(wq), w = 1-(ho/(b*b - b)), 同理画出区分50张原图。画出Figure.6
重采样：
Task1：计算一张重采样图片(resize 1.5)的p-map图
Task2：获得p-map图的频域图并滤波增强对比度
Task3：获得图像的特征值(频域能量突变值)
Task4：完成对100张图片及其1.5倍缩放图像的检测分类并完成ROC曲线
2、实验设备（环境）
Python，matlab，Opencv-python，numpy，matplotlib
3、实验步骤
中值滤波：这个实验是跟着曾辉老师做的。
（1）h0特征散点图，使用matlab编写如下代码
（2）rho=h0/h1特征散点图，使用matlab编写如下代码
（3）加权rho 特征散点图，使用matlab编写如下代码
重采样：
（4）计算一张重采样图片(resize 1.5)的p-map图，编写代码如下。
（5）获得p-map图的频域图并滤波增强对比度，编写代码如下。
（6）获得图像的特征值(频域能量突变值)，编写代码如下。
（7）先将100张图片1.5倍缩放之后保存起来。
（8）完成对100张图片及其1.5倍缩放图像的检测分类并完成ROC曲线，编写代码如下。
4、实验结果
中值滤波：这个实验是跟着曾辉老师做的。
（1）h0特征散点图
（2）rho=h0/h1特征散点图
（3）加权rho 特征散点图
（4）计算一张重采样图片(resize 1.5)的p-map图
这是1.5倍的原图
重采样图片(resize 1.5)的p-map图
（5）获得p-map图的频域图并滤波增强对比度
频域图
增强之后
（6）获得100张图片及其1.5倍缩放图像的特征值(频域能量突变值)
200张图
（7）完成对100张图片及其1.5倍缩放图像的检测分类并完成ROC曲线
5、讨论和分析（分析实验原理，为什么会得到这样的结果以及在试验中应注意的问题）
实验原理：中值滤波是一种非线性数字滤波器技术，经常用于去除图像或者其它信号中的噪声。这个设计思想就是检查输入信号中的采样并判断它是否代表了信号，使用奇数个采样组成的观察窗实现这项功能。观察窗口中的数值进行排序，位于观察窗中间的中值作为输出。然后，丢弃最早的值，取得新的采样，重复上面的计算过程。中值滤波是图像处理中的一个常用步骤，它对于斑点噪声和椒盐噪声来说尤其有用。保存边缘的特性使它在不希望出现边缘模糊的场合也很有用。
重采样：图像重采样包含两种情形，一种是下采样，把图像变小；另一种是上采样，把图像变大。重采样过程本质上是图像恢复过程,它用输入的离散数字图像重建代表原始图像二维连续函数,再按新的像元间距和像元位置进行采样。其数学过程是根据重建的连续函数(曲面),用周围若二像元点的值估计或内插出新采样点的值,相当于用采样函数与输入图像作二维卷积运算。
原始图像，经过中值滤波处理，图像更加趋于平滑，趋于一个平稳的值，去除噪声的同时保护信号的边缘。所以看出图像经过不断处理，最终趋于稳定，接近于一条直线。
图片放大后计算特征值和原图计算特征值画出ROC曲线，当把原始图像的特征值从大到小排序，不断和放大之后的特征值进行比较。发现，随着取的原始图像的特征值的不断减小，放大图像的特征值比其大的数量也越来越大，直到找到一个值可以使这两类图片完全区分开。
实验中应注意的问题：整个实验自己写起来代码还是有些困难的，在学长和老师帮助下解决一个一个小问题，最终成功实现。需要注意其中很多小问题，比如分号，小写L和1的区别等等。复习了Plot画图的过程，不仅可以画直线，也可以用点或者*来代替。学会如何批处理文件并将其放大或者缩小后重新保存。学会了如何借助python使用plt.plot()函数画出图像。
代码
chong_cai_yang
import numpy as np
import matplotlib.image as mpimg
import matplotlib.pyplot as plt
import os
import time
import cv2
def getcha(P,h):
P1 = P[:int(np.ceil(h/2)),int(np.ceil(h/2)):]
#公式23
C = np.array([0]*int(np.ceil(h/2))).astype(np.float32)
x=(np.sum(abs(P1)**2))**(-1)
for j in range(0,int(np.ceil(h/2))):
C[j]=np.sum(abs(P1[int(np.ceil(h/2))-1-j:-1,0:1+j])**2)*x
#进行分类
C = abs(C)
C = C / np.max(C)
cha=C[:-1]-C[1:]
return cha
start =time.perf_counter()
#读取文件
dirname=os.listdir("../data/resamping/image")
#设置初值
K=1
t=2
a=np.array([[-0.25,0.5,-0.25],[0.5,0,0.5],[-0.25,0.5,-0.25]])
#计算图像的特征值并进行分类
for name in dirname:
# print('%.2f%%' % (star/len(dirname)*100))
picdir = r'../data/resamping/image/'+name
img = mpimg.imread(picdir)
#img = img[:,:,0] ,如果img读取到三通道图，取一个通道来进行计算
h, w = img.shape
h2 = int(np.ceil(h / 2))
#对img矩阵四周添0以便于计算
padimg= np.pad(img,((1, 1),(1, 1)),'constant',constant_values=(0,0)).astype(np.float32)
#计算e
e=padimg-cv2.filter2D(padimg,-1,a)
e=e[1:-1,1:-1]
e = np.exp(-e ** 2)
#计算p
p=np.exp(-(abs(e)**2))
#对p进行离散傅里叶变换
V=np.fft.fft2(p)
V=np.fft.fftshift(V)
#对V进行滤波
mask = np.ones(img.shape,np.uint8)
win2 = 40/256*h2
mask[int(h/2-win2):int(h/2+win2),int(w/2-win2):int(w/2+win2)]=0
P0 = abs(V)*mask
P2 = (P0 / np.max(np.abs(P0))) ** 4 * np.max(np.abs(P0))
#计算特征值cha
cha = getcha(P2,h)
if abs(cha).max()>0.5:
print(abs(cha).max())
print("yes {}".format(name))
else:
print(abs(cha).max())
print("no {}".format(name))
end = time.perf_counter()
print('Running time: %s Seconds'%(end-start))roc
import matplotlib.pyplot as plt
import matplotlib
import numpy as np
import cv2
import os
import matplotlib.image as mping
import shutil
from PIL import Image
from scipy import ndimage
def getcha(P, h):
P0 = P[:int(np.ceil(h / 2)), int(np.ceil(h / 2)):]
l, m = np.where(P0 == np.max(P0))
m = m[0] + 1
if m < h / 2:
P1 = P0[:int(np.ceil(h / 2)), :]
else:
P1 = P0[:, int(np.ceil(h / 2)):]
C = np.array([0] * int(np.ceil(h / 2))).astype(np.float32)
x = (np.sum(abs(P1) ** 2)) ** (-1)
for j in range(0, int(np.ceil(h / 2))):
C[j] = np.sum(abs(P1[int(np.ceil(h / 2)) - 1 - j:-1, 0:1 + j]) ** 2) * x
C = abs(C)
C = C / np.max(C)
cha = C[:-1] - C[1:]
return cha
a = np.array([[-0.25, 0.5, -0.25], [0.5, 0, 0.5]
, [-0.25, 0.5, -0.25]])
r = 1
t = 2
K = 1
o = 1
dst = mping.imread(r"F:\1_need_delete\media_3_experiment\img\ucid00001.tif")
#dst=mping.imread("111.tif")
dst = dst[:, :,0]
h, w = np.shape(dst)
h2 = int(np.ceil(h / 2))
padimg = np.pad(dst, ((1, 1), (1, 1)), 'constant', constant_values=(0, 0)).astype(np.float32)
e = padimg - cv2.filter2D(padimg, -1, a)
e1 = padimg - cv2.filter2D(padimg, -1, a, borderType=cv2.BORDER_CONSTANT)
e = e[1:-1, 1:-1]
e = np.exp(-e ** 2)
p = r * np.exp(-(abs(e) ** 2) / o)
#cv2.imwrite("4444.tif", p)
V = np.fft.fft2(p)
V = np.fft.fftshift(V)
mask = np.ones(dst.shape, np.uint8)
win2 = 40 / 256 * h2
mask[int(h / 2 - win2):int(h / 2 + win2), int(w / 2 - win2):int(w / 2 + win2)] = 0
P0 = abs(V) * mask
P2 = ((P0 / np.max(np.abs(P0))) ** 4) * (np.max(np.abs(P0)))
print(np.where(P2 == P2.max()))
plt.imshow(P2)
plt.show()seniors
# -*- coding:utf-8 -*-
import numpy as np
import cv2
import matplotlib.pyplot as plt
a=np.array([[-0.25,0.5,-0.25],[0.5,0,0.5],[-0.25,0.5,-0.25]])
k=1
r=1
t=2
o=1
def show(img):
img = cv2.cvtColor(img, cv2.COLOR_BGR2RGB)
plt.imshow(img)
plt.show()
def getPmap(img2):
pading=np.pad(img2,((1,1),(1,1)),'constant',constant_values=(0,0)).astype(np.float32)
##cv2.filter2D(img2,-1,a,borderType=cv2.BORDER_CONSTANT)
e=pading-cv2.filter2D(pading,-1,a)
e=e[1:-1,1:-1]
e=np.exp(-e**2)
p=r*np.exp(-(abs(e)**2)/o)
return p
def FFT(p,img):
mask=np.zeros(img.shape,)
P1=np.fft.fft2(p)
P1=np.fft.fftshift(P1)
P0=abs(P1)*mask
if __name__=="__main__":
for i in range(1,101):
n=str(i)
src = "img/ucid" + n.zfill(5) + ".tif"
img=cv2.imread(src,0)
(h,w)=img.shape
img2=cv2.resize(img,(int(1.5*w),int(1.5*h)))
show(img2)
result=getPmap(img2)
show(result)
getPmap(img)
break
