其他
rgee学习笔记之Image（二）
rgee学习笔记之Image（二）
前面rgee学习笔记之Image（一）介绍了rgee中Image的展示，以及元数据等基本影像信息的获取。下面继续学习Image的用法。
Mathematical operations数值运算
这部分代码主要是进行rgee的波段运算。
加载Landsat7数据
# Load two 5-year Landsat 7 composites.
landsat1999 <- ee$Image("LANDSAT/LE7_TOA_5YEAR/1999_2003")
landsat2008 <- ee$Image("LANDSAT/LE7_TOA_5YEAR/2008_2012")
使用复杂的方法计算NDVI
# Compute NDVI the hard way.
ndvi1999 <- landsat1999$select("B4")$subtract(landsat1999$select("B3"))$divide(landsat1999$select("B4")$add(landsat1999$select("B3")))
使用GEE内置函数计算NDVI normalizedDifference是GEE中的内置函数，查阅参考文献2中的文档可知：计算公式：(first − second) / (first + second) 输入：波段名 输出：Image
# Compute NDVI the easy way.
ndvi2008 <- landsat2008$normalizedDifference(c("B4", "B3"))
Image.subtract(image2)GEE内置函数，相减
# Compute the multi-band difference image.
diff <- landsat2008$subtract(landsat1999)
Image.pow(image2)GEE内置函数，幂运算
# Compute the squared difference in each band.
squaredDifference <- diff$pow(2)
计算结果的显示
Map$setZoom(zoom = 3)
Map$addLayer(
eeObject = diff,
visParams = list(bands = c("B4", "B3", "B2"), min = -32, max = 32),
name = "difference"
) +
Map$addLayer(
eeObject = squaredDifference,
visParams = list(bands = c("B4", "B3", "B2"), max = 1000),
name = "squared diff."
)
Relational, conditional and Boolean operations关系、条件和布尔运算
Image.lt(image2)当image1 < image2时，输出1输入是image 输出是布尔值 Image.and(image2)当image1和image2都不为0时，返回1输入是image 输出是布尔值
# Load a Landsat 8 image.
image <- ee$Image("LANDSAT/LC08/C01/T1_TOA/LC08_044034_20140318")
# Create NDVI and NDWI spectral indices.
ndvi <- image$normalizedDifference(c("B5", "B4"))
ndwi <- image$normalizedDifference(c("B3", "B5"))
# Create a binary layer using logical operations.
bare <- ndvi$lt(0.2)$And(ndwi$lt(0))
# Mask and display the binary layer.
Map$setCenter(lon = -122.3578, lat = 37.7726)
Map$setZoom(zoom = 12)
Map$addLayer(
eeObject = bare$updateMask(bare),
visParams = list(min = 0, max = 20),
name = "bare"
)
Convolutions卷积
ee.Kernel.square
`ee.Kernel.square(radius, units, normalize, magnitude) radius,浮点型，核半径 units，字符串，默认"pixels"，可以是"pixels"或者"meters"，指定核半径单位，如果核半径单位是米，当缩放级别（zoom-level）变化的时候会变化 normalize，布尔值，默认："true"，标准化核值和为1（Normalize the kernel values to sum to 1） magnitude，浮点型，默认1，确定每个值数量级别（Scale each value by this amount）
ee.Image.convolve
Image.convolve(kernel)利用给出的核（Kernel）进行卷积运算图像卷积操作的目的是利用像素点和其邻域像素之前的空间关系，通过加权求和的操作，实现模糊（blurring），锐化（sharpening），边缘检测（edge detection）等功能。 卷积的概念和常见运算可以查阅参考文献5
ee.Kernel.laplacian
`ee.Kernel.laplacian8(magnitude, normalize)`` 生成一个3×3的Laplacian-8边缘检测核
# Load and display an image.
image <- ee$Image("LANDSAT/LC08/C01/T1_TOA/LC08_044034_20140318")
# Define a boxcar or low-pass kernel.
# boxcar = ee$Kernel$square(list(
# radius = 7, units = 'pixels', 'normalize' = T
# ))
boxcar <- ee$Kernel$square(7, "pixels", T)
# Smooth the image by convolving with the boxcar kernel.
smooth <- image$convolve(boxcar)
# Define a Laplacian, or edge-detection kernel.
laplacian <- ee$Kernel$laplacian8(1, F)
# Apply the edge-detection kernel.
edgy <- image$convolve(laplacian)
Map$setCenter(lon = -121.9785, lat = 37.8694)
Map$setZoom(zoom = 11)
Map$addLayer(
eeObject = image,
visParams = list(bands = c("B5", "B4", "B3"), max = 0.5),
name = "input image"
) +
Map$addLayer(
eeObject = smooth,
visParams = list(bands = c("B5", "B4", "B3"), max = 0.5),
name = "smoothed"
) +
Map$addLayer(
eeObject = edgy,
visParams = list(bands = c("B5", "B4", "B3"), max = 0.5),
name = "edges"
)
Morphological Operations形态操作
ee.Kernel.circle
`ee.Kernel.circle(radius, units, normalize, magnitude)`` 生成一个圆形的布尔核
ee.Image.focal_min
使用一个定义的核，应用一个形态学的滤波 `Image.focal_min(radius, kernelType, units, iterations, kernel)`` radius, 浮点型，默认1.5 kernelType，字符串，默认“circle”，可选的包括：'circle', 'square', 'cross', 'plus', octagon' and 'diamond' units,字符串，默认"pixels"，可选"meters"或"pixels" iterations，整型，默认1，应用核的次数 kernel，核类型，默认：空，一个自定义的核，如果定义了核，kernelType和radius可以忽略 类似的函数有： ee.Image.focal_medianee.Image.focal_meanee.Image.focal_max
# Load a Landsat 8 image, select the NIR band, threshold, display.
# 加载Landsat8影像，选择近红外波段，大于0.2的部分
image <- ee$Image("LANDSAT/LC08/C01/T1_TOA/LC08_044034_20140318")$select(4)$gt(0.2)
# Define a kernel.
kernel <- ee$Kernel$circle(radius = 1)
# Perform an erosion followed by a dilation, display.
opened <- image$focal_min(kernel = kernel, iterations = 2)$focal_max(kernel = kernel, iterations = 2)
Map$setCenter(lon = -122.1899, lat = 37.5010)
Map$setZoom(zoom = 13)
Map$addLayer(
eeObject = image,
visParams = list(min = 0, max = 1),
name = "NIR threshold"
) +
Map$addLayer(
eeObject = opened,
visParams = list(min = 0, max = 1),
name = "opened"
)
Gradients梯度
Image.gradient()计算x,y梯度梯度一般用于增强图像对比度 梯度是一个向量 梯度的方向是函数变化最快的方向 图像边缘处变化大，灰度值变化大，梯度大 图像平滑处变化小，灰度值变化小，梯度小
# Load a Landsat 8 image and select the panchromatic band.
image <- ee$Image("LANDSAT/LC08/C01/T1/LC08_044034_20140318")$select("B8")
# Compute the image gradient in the X and Y directions.
xyGrad <- image$gradient()
# Compute the magnitude of the gradient.
gradient <- xyGrad$select("x")$pow(2)$add(xyGrad$select("y")$pow(2))$sqrt()
# Compute the direction of the gradient.
direction <- xyGrad$select("y")$atan2(xyGrad$select("x"))
# Display the results.
Map$setCenter(lon = -122.054, lat = 37.7295)
Map$setZoom(zoom = 10)
Map$addLayer(
eeObject = direction,
visParams = list(min = -2, max = 2),
name = "direction"
) +
Map$addLayer(
eeObject = gradient,
visParams = list(min = -7, max = 7),
name = "opened"
)
参考文献
https://csaybar.github.io/rgee-examples/ https://developers.google.com/earth-engine/api_docs/#eeimagenormalizeddifference https://developers.google.com/earth-engine/api_docs/#eeimagesubtract https://developers.google.com/earth-engine/api_docs/#eekernelsquare https://blog.csdn.net/u012005313/article/details/84068337 https://developers.google.com/earth-engine/api_docs/#eeimagefocal_min https://zhuanlan.zhihu.com/p/113397988