landsatL2SP数据还用辐射定标标定吗？

来源：蜘蛛抓取(WebSpider) 时间：2023-03-21 10:38 标签： 辐射定标

加载包；
library(pacman)
p_load(raster,sp,rgdal,rgeos,tidyverse)
加载影像，并以真彩色显示；
img <- stack('Anhui/LC08_L2SP_121038_20201024_20201105_02_T1/hfLandsat8.tif')
plotRGB(img,r=3,g=2,b=1,axes=TRUE,stretch='lin',main="Landsat True Colour Composite")
nlayers(img)#查看图层数
res(img)#查看空间分辨率
img
# class
: RasterStack
# dimensions : 5884, 4034, 23736056, 6
(nrow, ncol, ncell, nlayers)
# resolution : 30, 30
(x, y)
# extent
: 466845, 587865, 3424335, 3600855
(xmin, xmax, ymin, ymax)
# crs
: +proj=utm +zone=50 +datum=WGS84 +units=m +no_defs
# names
: hfLandsat8.1, hfLandsat8.2, hfLandsat8.3, hfLandsat8.4, hfLandsat8.5, hfLandsat8.6
# min values :
-31788,
-32711,
-32103,
-32762,
-32745,
-32748
# max values :
32502,
32294,
32749,
32746,
32749,
32741
把raster转为data frame
img.df <- as.data.frame(img)#把raster转为dataframe
str(img.df)#查看数据结构
# 'data.frame': 23736056 obs. of
6 variables:
# $ hfLandsat8.1: num
NA NA NA NA NA NA NA NA NA NA ...
# $ hfLandsat8.2: num
NA NA NA NA NA NA NA NA NA NA ...
# $ hfLandsat8.3: num
NA NA NA NA NA NA NA NA NA NA ...
# $ hfLandsat8.4: num
NA NA NA NA NA NA NA NA NA NA ...
# $ hfLandsat8.5: num
NA NA NA NA NA NA NA NA NA NA ...
# $ hfLandsat8.6: num
NA NA NA NA NA NA NA NA NA NA ...
colnames(img.df)<-c('Band1','Band2','Band3','Band4','Band5','Band6')#重命名dataframe的列
head(img.df)#查看前几行
# Band1 Band2 Band3 Band4 Band5 Band6
# 1
NA
NA
NA
NA
NA
NA
# 2
NA
NA
NA
NA
NA
NA
# 3
NA
NA
NA
NA
NA
NA
# 4
NA
NA
NA
NA
NA
NA
# 5
NA
NA
NA
NA
NA
NA
# 6
NA
NA
NA
NA
NA
NA
提取经纬度到dataframe
SpatialPoints(img,CRS("+init=epsg:32650"))%>%
spTransform(., CRS("+init=epsg:4326"))->lonlat1
head(lonlat1)
lonlat1.df <-as.data.frame(lonlat1,na.rm=TRUE)
head(lonlat1.df)
colnames(lonlat1.df)<-c('Lon','Lat')
head(lonlat1.df)
df<-cbind(img.df,lonlat1.df)
head(df)
# Band1 Band2 Band3 Band4 Band5 Band6
Lon
Lat
# 1
NA
NA
NA
NA
NA
NA 116.6470 32.54444
# 2
NA
NA
NA
NA
NA
NA 116.6474 32.54444
# 3
NA
NA
NA
NA
NA
NA 116.6477 32.54444
# 4
NA
NA
NA
NA
NA
NA 116.6480 32.54444
# 5
NA
NA
NA
NA
NA
NA 116.6483 32.54444
# 6
NA
NA
NA
NA
NA
NA 116.6486 32.54444
计算NDVI
NDVICalc <- function(x){
ndvi<-(x[[4]]-x[[3]])/(x[[3]]+x[[4]])
return(ndvi)
}
NDVI <- NDVICalc(img)
plot(NDVI)
很明显上面计算的NDVI中有异常值，我们NDVI应该是在-1至1之间的；怎么办呢？我也不知道。。。Google下吧！果然找到了解决方法,先看下NDVI的直方图，和NDVI信息（values : -31382, 109.7899 (min, max)）
hist(NDVI,main="NDVI_hist",xlab="NDVI_VALUE",ylab="Frequency",col="wheat1")
NDVI
# class
: RasterLayer
# dimensions : 5884, 4034, 23736056
(nrow, ncol, ncell)
# resolution : 30, 30
(x, y)
# extent
: 466845, 587865, 3424335, 3600855
(xmin, xmax, ymin, ymax)
# crs
: +proj=utm +zone=50 +datum=WGS84 +units=m +no_defs
# source
: memory
# names
: layer
# values
: -31382, 109.7899
(min, max)
利用阈值法把-1,1范围外的赋为NA，然后再看下NDVI的原信息这回，values : -0.8023534, 0.8727696 (min, max)，搞定。
NDVI[NDVI>1|NDVI<=-1]<-NA
NDVI
hist(NDVI)
# class
: RasterLayer
# dimensions : 5884, 4034, 23736056
(nrow, ncol, ncell)
# resolution : 30, 30
(x, y)
# extent
: 466845, 587865, 3424335, 3600855
(xmin, xmax, ymin, ymax)
# crs
: +proj=utm +zone=50 +datum=WGS84 +units=m +no_defs
# source
: memory
# names
: layer
# values
: -0.8023534, 0.8727696
(min, max)
plot(NDVI)
把NDVI存储为dataframe，并统计每列
NDVI.df <- as.data.frame(NDVI)
str(NDVI.df)
colnames(NDVI.df)<-"NDVI"
head(NDVI.df)
allvars <- cbind(df,NDVI.df)
head(allvars)
allvars
allvars <-na.exclude(allvars)
summary(allvars)
# Band1
Band2
Band3
Band4
Band5
Band6
# Min.
:-31788
Min.
:-32711
Min.
:-32103
Min.
:-29310
Min.
:-32745
Min.
:-32748
# 1st Qu.:
8786
1st Qu.:
9546
1st Qu.:
9382
1st Qu.: 13686
1st Qu.: 12085
1st Qu.:
9993
# Median :
9167
Median : 10040
Median : 10089
Median : 15298
Median : 13600
Median : 11164
# Mean
:
9289
Mean
: 10169
Mean
: 10216
Mean
: 14641
Mean
: 13252
Mean
: 11190
# 3rd Qu.:
9660
3rd Qu.: 10667
3rd Qu.: 10949
3rd Qu.: 16374
3rd Qu.: 14905
3rd Qu.: 12459
# Max.
: 32502
Max.
: 32294
Max.
: 29784
Max.
: 32746
Max.
: 32749
Max.
: 32741
# Lon
Lat
NDVI
# Min.
:116.6
Min.
:30.95
Min.
:-0.8024
# 1st Qu.:117.1
1st Qu.:31.53
1st Qu.: 0.1255
# Median :117.3
Median :31.77
Median : 0.1897
# Mean
:117.3
Mean
:31.76
Mean
: 0.1693
# 3rd Qu.:117.5
3rd Qu.:32.01
3rd Qu.: 0.2488
# Max.
:117.9
Max.
:32.54
Max.
: 0.8728

我要回帖

更多关于 辐射定标 的文章

 

随机推荐