1. 程式人生 > >python環境下的遙感影象dct水印提取

python環境下的遙感影象dct水印提取

插入在前面的文章 本文主要是提取演算法,詳細註釋看插入演算法import gdalimport numpyfrom skimage import data,img_as_float,img_as_ubyteimport cv2filename = "f:/bs/band_11_255_0.0075.tif"#插入水印後的圖片dataset = gdal.Open(filename)im_width = dataset.RasterXSize #獲取長im_height = dataset.RasterYSize #獲取寬im_bands = dataset.RasterCount #獲取通道數im_data = dataset.ReadAsArray(0,0,im_width,im_height)BLOCK = numpy.zeros([32,32])M = 3808N = 238K = 16I = numpy.zeros([M,M])I = im_data[2000:,2000:]filename = "f:/bs/band_1.tif"dataset = gdal.Open(filename)im_width = dataset.RasterXSize #獲取長im_height = dataset.RasterYSize #獲取寬im_bands = dataset.RasterCount #獲取通道數P = dataset.ReadAsArray(0,0,im_width,im_height)P = P[2000:,2000:]W = numpy.zeros([N,N])for p in range(238): for q in range (238): x = p*K y = q*K BLOCK1 = P[x:x+K,y:y+K] BLOCK2 = I[x:x+K,y:y+K] BLOCK1 = img_as_float(BLOCK1) #轉換unit8為float BLOCK2 = img_as_float(BLOCK2) BLOCK1 = cv2.dct(BLOCK1) BLOCK2 = cv2.dct(BLOCK2) a = BLOCK1[0][0]/BLOCK2[0][0]-1 if a<0: W[p][q]=1 else : W[p][q]=-1W = img_as_ubyte(W) #轉換float為unit8format = "GTiff"driver = gdal.GetDriverByName(format)tods = driver.Create("f:/bs/water_wh_1_255_0.0075.bmp",238,238,1,options=["INTERLEAVE = PIXEL"])tods.WriteRaster(0,0,238,238,W.tostring(),238,238,band_list=[1])原圖:插入水印後:水印:提取出的水印: