超碰日本,亚洲在线成人视频,一级A片电影A片录像,91高潮久久久久久久,www久久,色婷婷在线视频网站,黄网站免费看欧美,日韩99在线观看

點(diǎn)擊上方“小白學(xué)視覺”，選擇加"星標(biāo)"或“置頂”

重磅干貨，第一時(shí)間送達(dá)

介紹

海岸是一個(gè)動(dòng)態(tài)系統(tǒng)，其中因侵蝕現(xiàn)象導(dǎo)致的海岸線后退、或是由眾多因素如氣象，地質(zhì)，生物和人類活動(dòng)所導(dǎo)致線前進(jìn)的是常見現(xiàn)象。

在海洋磨損作用大于沉積物的情況下，有明顯的海岸侵蝕，我們稱之為地球表面的崩解和破壞。

資料來源：弗林德斯大學(xué)（CC0）

本文的目標(biāo)

在本文中，我們將對(duì)Landsat 8平臺(tái)上的OLI（陸地成像儀）傳感器獲取的衛(wèi)星圖像使用Canny Edge Detection算法。

通過這種方法，我們將能夠可視化的估計(jì)特定歐洲地區(qū)遭受強(qiáng)腐蝕作用的海岸線隨時(shí)間的推移：霍德內(nèi)斯海岸。

一下是處理流程：

處理流程

在開始之前讓我們先介紹一下OLI數(shù)據(jù)...

0.關(guān)于Landsat OLI數(shù)據(jù)的簡要介紹

Landsat 8是一個(gè)軌道平臺(tái)，安裝在稱為OLI（陸地成像儀）的11波段多光譜傳感器上。

具體來說，在本文中，我們將僅使用分辨率為30米（即前7個(gè)）的波段。

美國地質(zhì)調(diào)查局陸地衛(wèi)星8號(hào)

該數(shù)據(jù)可以免費(fèi)下載，注冊后，獲得USGS：https://earthexplorer.usgs.gov/。

而且，通常我摸并不使用入射太陽光作為原始數(shù)據(jù)，而是使用反射率，即從地球表面反射的太陽光量[0-1]。

1.包導(dǎo)入

在各種常見的包，我們將使用rasterio處理圖像，利用OpenCV中的Canny?算法和Scikit-Learn分割圖像。

from glob import globimport numpy as np
import rasterioimport json, re, itertools, os
import matplotlib.pyplot as plt
import cv2 as cvfrom sklearn import preprocessingfrom sklearn.cluster import KMeans

2.數(shù)據(jù)導(dǎo)入

讓我們定義一個(gè)變量，該變量告訴我們要保留的波段數(shù)以及在JSON中輸入的輔助數(shù)據(jù)：

N_OPTICS_BANDS = 7
with open("bands.json","r") as bandsJson:    bandsCharacteristics = json.load(bandsJson)

這個(gè)Json是Landsat OLI成像儀的信息集合。類似于一種說明手冊：

# bands.json
[{'id': '1', 'name': 'Coastal aerosol', 'span': '0.43-0.45', 'resolution': '30'}, {'id': '2', 'name': 'Blue', 'span': '0.45-0.51', 'resolution': '30'}, {'id': '3', 'name': 'Green', 'span': '0.53-0.59', 'resolution': '30'}, {'id': '4', 'name': 'Red', 'span': '0.64-0.67', 'resolution': '30'}, {'id': '5', 'name': 'NIR', 'span': '0.85-0.88', 'resolution': '30'}, {'id': '6', 'name': 'SWIR 1', 'span': '1.57-1.65', 'resolution': '30'}, {'id': '7', 'name': 'SWIR 2', 'span': '2.11-2.29', 'resolution': '30'}, {'id': '8', 'name': 'Panchromatic', 'span': '0.50-0.68', 'resolution': '15'}, {'id': '9', 'name': 'Cirrus', 'span': '1-36-1.38', 'resolution': '30'}, {'id': '10', 'name': 'TIRS 1', 'span': '10.6-11.9', 'resolution': '100'}, {'id': '11', 'name': 'TIRS 2', 'span': '11.50-12.51', 'resolution': '100'}]

bands.json文件包含有關(guān)我們將要使用的頻段的所有有用信息。

注意，我們將僅使用分辨率為30 m的頻段，因此僅使用前7個(gè)頻段。如果您愿意使用較低的分辨率（100m），則也可以嵌入TIRS 1和TIRS 2頻段。

正如上面幾行已經(jīng)提到的那樣，我們將使用從Landsat-8 OLI上獲取兩組不同的數(shù)據(jù)：

? 2014/02/01

? 2019/07/25

為了簡化兩次采集的所需操作，我們將定義一個(gè)Acquisition（）類，其中將封裝所有必要的函數(shù)。

在執(zhí)行代碼期間，我們能夠執(zhí)行一些基礎(chǔ)支持性的功能，例如：

? 在指定路徑中搜索GeoTIFF；

? 加載采購；

? 購置登記?（調(diào)整）；

? 收購子集

class Acquisition:        def __init__(self, path, ext, nOpticsBands):        self.nOpticsBands = nOpticsBands        self._getGeoTIFFs(path, ext)        self.images = self._loadAcquisition()            def _getGeoTIFFs(self, path, ext):        # It searches for GeoTIFF files within the folder.        print("Searching for '%s' files in %s" % (ext, path))       self.fileList = glob(os.path.join(path,"*."+ext))        self.opticsFileList = [    [list(filter(re.compile(r"band%s\."%a).search, self.fileList))[0] for a in range(1,self.nOpticsBands+1)]        print("Found %d 'tif' files" % len(self.opticsFileList))            def _loadAcquisition(self):        # It finally reads and loads selected images into arrays.        print("Loading images")        self.loads = [rasterio.open(bandPath) for bandPath in self.opticsFileList]        images = [load.read()[0] for load in self.loads]        print("Done")        return images        def subsetImages(self, w1, w2, h1, h2, leftBound):        # This function subsets images according the defined sizes.        print("Subsetting images (%s:%s, %s:%s)" % (w1, w2, h1, h2))        cols = (self.loads[0].bounds.left - leftBound)/30        registered = [np.insert(band,np.repeat(0,cols),0,axis=1) for band in self.images]        subset = [band[w1:w2,h1:h2] for band in registered]        print("Done")        return subset

好的，讓我們現(xiàn)在開始啟動(dòng)整個(gè)代碼：

DATES = ["2014-02-01", "2019-07-25"]acquisitionsObjects = []
for date in DATES:    singleAcquisitionObject = Acquisition("Data/"+date, "tif", N_OPTICS_BANDS)    acquisitionsObjects.append( singleAcquisitionObject )

運(yùn)行結(jié)果如下：

Searching for 'tif' files in Data/2014-02-01
Found 7 'tif' files
Loading images
Done
Searching for 'tif' files in Data/2019-07-25
Found 7 'tif' files
Loading images
Done

現(xiàn)在我們已加載了14張OLI圖像（在7個(gè)波段中各采集2個(gè)）。

2.1 子集多光譜立方體

在這個(gè)階段中，先對(duì)兩個(gè)多光譜立方體進(jìn)行“對(duì)齊”（或正式注冊），再切出不感興趣的部分。

我們可以使用ImageImages（）函數(shù)“剪切”不需要的數(shù)據(jù)。

因此，我們定義AOI（感興趣的區(qū)域），并使用Acquisition（）類中的subsetImages（）函數(shù)進(jìn)行設(shè)置：

W1, W2 = 950, 2300H1, H2 = 4500, 5300
subAcquisitions = [acquisition.subsetImages(W1, W2, H1, H2, 552285.0) for acquisition in acquisitionsObjects].

完成！

3.數(shù)據(jù)探索

3.1可視化多光譜立方體

讓我們嘗試查看2019/07/25收購的所有范圍。出于純粹的美學(xué)原因，在繪制圖像之前，讓我們使用

StandardScaler（）對(duì)圖像進(jìn)行標(biāo)準(zhǔn)化。

axs = range(N_OPTICS_BANDS)fig, axs = plt.subplots(2, 4, figsize=(15,12))axs = list(itertools.chain.from_iterable(axs))
for b in range(N_OPTICS_BANDS):    id_ = bandsCharacteristics[b]["id"]    name_ = bandsCharacteristics[b]["name"]    span_ = bandsCharacteristics[b]["span"]    resolution_ = bandsCharacteristics[b]["resolution"]    title = "%s - %s\n%s (%s m)" % (id_, name_, span_, resolution_)    axs[b].imshow(preprocessing.StandardScaler().fit_transform(subAcquisitions[1][b]), cmap="Greys_r")    axs[b].set_title(title); axs[b].set_xticklabels([]); axs[b].set_yticklabels([])
plt.axis("off"); plt.tight_layout(w_pad=-10); plt.show()

以下是運(yùn)行結(jié)果。

這些圖中，有些波段比其他波段更亮。這很正常。

3.2可視化復(fù)合RGB中的多光譜立方體

現(xiàn)在，讓我們嘗試可視化使用波段4（紅色），3（綠色）和2（藍(lán)色）獲得的RGB復(fù)合圖像中的兩次采集。

定義BIAS和GAIN?僅是為了獲得更好的效果。

BIAS = 1.5GAIN = [2.3,2.4,1.4]
r1 = (subAcquisitions[0][3] - subAcquisitions[0][3].min()) / (subAcquisitions[0][3].max()-subAcquisitions[0][3].min()) * GAIN[0] * BIASg1 = (subAcquisitions[0][2] - subAcquisitions[0][2].min()) / (subAcquisitions[0][2].max()-subAcquisitions[0][2].min()) * GAIN[1] * BIASb1 = (subAcquisitions[0][1] - subAcquisitions[0][1].min()) / (subAcquisitions[0][1].max()-subAcquisitions[0][1].min()) * GAIN[2] * BIAS
r2 = (subAcquisitions[1][3] - subAcquisitions[1][3].min()) / (subAcquisitions[1][3].max()-subAcquisitions[1][3].min()) * GAIN[0] * BIASg2 = (subAcquisitions[1][2] - subAcquisitions[1][2].min()) / (subAcquisitions[1][2].max()-subAcquisitions[1][2].min()) * GAIN[1] * BIASb2 = (subAcquisitions[1][1] - subAcquisitions[1][1].min()) / (subAcquisitions[1][1].max()-subAcquisitions[1][1].min()) * GAIN[2] * BIAS
rgbImage1, rgbImage2 = np.zeros((W2-W1,H2-H1,3)), np.zeros((W2-W1,H2-H1,3))rgbImage1[:,:,0], rgbImage2[:,:,0] = r1, r2rgbImage1[:,:,1], rgbImage2[:,:,1] = g1, g2rgbImage1[:,:,2], rgbImage2[:,:,2] = b1, b2
fig, (ax1,ax2) = plt.subplots(1,2,figsize=(16,12))ax1.imshow(rgbImage1); ax2.imshow(rgbImage2)ax1.set_title("RGB\n(Bands 4-3-2)\n2014-02-01"); ax2.set_title("RGB\n(Bands 4-3-2)\n2019-07-25")plt.show()

結(jié)果如下圖所示！有趣的是，這兩次獲取的反射率完全的不同。

好的，繼續(xù)進(jìn)行海岸線檢測。

4.自動(dòng)化海岸線檢測

在本段中，我們將使用Canny的算法執(zhí)行邊緣檢測。

在進(jìn)行實(shí)際檢測之前，有必要準(zhǔn)備數(shù)據(jù)集，嘗試通過聚類算法對(duì)數(shù)據(jù)集進(jìn)行分割以區(qū)分海洋和陸地。

4.1數(shù)據(jù)準(zhǔn)備

在此階段，我們將重塑兩個(gè)多光譜立方體以進(jìn)行聚類操作。

4.2用K均值進(jìn)行圖像分割

我們通過k均值對(duì)這兩次采集進(jìn)行細(xì)分（使用自己喜歡的模型即可）。

4.3細(xì)分結(jié)果

這是確定的代表新興土地和水體的兩個(gè)集群。

4.4Canny邊緣檢測算法

Canny的傳統(tǒng)鍵技術(shù)分為以下幾個(gè)階段：

1. 高斯濾波器通過卷積降低噪聲；

2. 四個(gè)方向（水平，垂直和2個(gè)傾斜）的圖像梯度計(jì)算；

3. 梯度局部最大值的提取;

4. 帶有滯后的閾值，用于邊緣提取。

讓我們開始，將聚類結(jié)果轉(zhuǎn)換為圖像，然后通過具有15x15內(nèi)核的高斯濾波器降低噪聲：

clusteredImages = [clusterLabels.reshape(subAcquisitions[0][0].shape).astype("uint8") for clusterLabels in clusters]blurredImages = [cv.GaussianBlur(clusteredImage, (15,15), 0) for clusteredImage in clusteredImages]
fig, (ax1, ax2) = plt.subplots(1,2,figsize=(16,13))
ax1.imshow(blurredImages[0])ax1.set_title("2014-02-01\nGaussian Blurred Image")ax2.imshow(blurredImages[1])ax2.set_title("2019-07-25\nGaussian Blurred Image")
plt.show()

在圖像稍微模糊之后，我們可以使用OpenCV? Canny（）模塊：

rawEdges = [cv.Canny(blurredImage, 2, 5).astype("float").reshape(clusteredImages[0].shape) for blurredImage in blurredImages]
edges = []for edge in rawEdges:    edge[edge == 0] = np.nan    edges.append(edge)

在單行代碼中，我們獲得了梯度，提取了局部最大值，然后對(duì)每次采集都應(yīng)用了帶有滯后的閾值。

注意：我們可以使用不同參數(shù)Canny（）來探索處理結(jié)果。

4.5結(jié)果

plt.figure(figsize=(16,30))plt.imshow(rgbImage2)plt.imshow(edges[0], cmap = 'Set3_r')plt.imshow(edges[1], cmap = 'Set1')plt.title('CoastLine')plt.show()

以下是一些詳細(xì)信息：

5結(jié)論

從結(jié)果中可以看到，Canny的算法在其原始管道中運(yùn)行良好，但其性能通常取決于所涉及的數(shù)據(jù)。

實(shí)際上，所使用的聚類算法使我們能夠?qū)Χ喙庾V立方體進(jìn)行細(xì)分。并行使用多個(gè)聚類模型可以總體上改善結(jié)果。

交流群

歡迎加入公眾號(hào)讀者群一起和同行交流，目前有SLAM、三維視覺、傳感器、自動(dòng)駕駛、計(jì)算攝影、檢測、分割、識(shí)別、醫(yī)學(xué)影像、GAN、算法競賽等微信群（以后會(huì)逐漸細(xì)分），請(qǐng)掃描下面微信號(hào)加群，備注：”昵稱+學(xué)校/公司+研究方向“，例如：”張三?+?上海交大?+?視覺SLAM“。請(qǐng)按照格式備注，否則不予通過。添加成功后會(huì)根據(jù)研究方向邀請(qǐng)進(jìn)入相關(guān)微信群。請(qǐng)勿在群內(nèi)發(fā)送廣告，否則會(huì)請(qǐng)出群，謝謝理解~

基于OpenCV實(shí)現(xiàn)海岸線變化檢測