-
Notifications
You must be signed in to change notification settings - Fork 0
/
Copy pathR_code_remote sensing.r
214 lines (147 loc) · 8.86 KB
/
R_code_remote sensing.r
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
# R code for Remote Sensing
setwd ("C:/lab/")
library (raster) # install.packages ("raster")
install.packages("RStoolbox") #packages to make our analysis, the book that the Prof. wrote is about the algorith of this R function that he created
library (RStoolbox) #to be faster in the istallation of pkgs: install.packages(c("raster", "RStoolbox")
#load images from the lab folder
p224r63_2011 <- brick ("p224r63_2011_masked.grd")
plot(p224r63_2011)
#landsat resolution of 30m (each pixel)...B1 blu, B2 green, B3 red, B4 NIR ecc
# B1: blue
# B2: green
# B3: red
# B4: near infrared (nir)
# B5: medium infrared
# B6: thermal infrared
# B7: medium infrared
ls()
cl <- colorRampPalette(c('black','grey','light grey'))(100) #changing our scaled colors for graphs (gray)
plot(p224r63_2011, col=cl) #Exercise: plot the image with scale color ramp palette built
cllow <- colorRampPalette(c('black','grey','light grey'))(5) # grey scaled low amount of colours
plot(p224r63_2011, col=cllow)
names(p224r63_2011) # [1] "B1_sre" "B2_sre" "B3_sre" "B4_sre" "B5_sre" "B6_bt" "B7_sre"
clb <- colorRampPalette(c('dark blue','blue','light blue'))(100) ##changing our scaled colors for graphs (blue)
plot(p224r63_2011$B1_sre, col=clb) #symble $ links together collum with the data I want
clnir <- colorRampPalette(c('red','orange','yellow'))(100) # Exercise: plot the near infrared band with:
plot(p224r63_2011$B4_sre, col=clnir) #colorRampPalette color that go from red, to orange, until yellow
dev.off ()
# multiframe to plot more graphs together, we say at the console how many row and collums we want in the Plot Space: par(mfrow=c(2,2))
par(mfrow=c(2,2))
# B1 blue band, we want set colors intensity in blue
clb <- colorRampPalette(c('dark blue','blue','light blue'))(100)
plot(p224r63_2011$B1_sre, col=clb) # func. par, multivariateframe = col = c finally names the colors you would
# plot the immages with colors calls clb
#dev.off to close all the windows if I need
# Exercise: do the same for the B2_sre green band
clg <- colorRampPalette(c('dark green','green','light green'))(100)
plot(p224r63_2011$B1_sre, col=clg)
#B3_sre red band
clr <- colorRampPalette(c('dark red','red','pink'))(100)
plot(p224r63_2011$B1_sre, col=clr)
#B4_sre near infrared band
cln <- colorRampPalette(c('red','orange','yellow'))(100)
plot(p224r63_2011$B1_sre, col=cln)
#change the layout of our graph in 4 rows and 1 collum and to make it different and stretched, could be really usefull for data analysis
dev.off ()
par(mfrow=c(4,1))
# B1: blue
clb <- colorRampPalette(c('dark blue','blue','light blue'))(100)
plot(p224r63_2011$B1_sre, col=clb)
# B2 green
# Exercise: do the same for the green band B2_sre
clg <- colorRampPalette(c('dark green','green','light green'))(100)
plot(p224r63_2011$B2_sre, col=clg)
# B3 red
clr <- colorRampPalette(c('dark red','red','pink'))(100)
plot(p224r63_2011$B3_sre, col=clr)
# B4 NIR
cln <- colorRampPalette(c('red','orange','yellow'))(100)
plot(p224r63_2011$B4_sre, col=cln)
# natural colours
# 3 components: R G B
# 3 bands: R = red band, G = green band, B = blue band
# plotRGB(p224r63_2011,r=3,g=2,b=1)
# B1: blue - 1
# B2: green - 2
# B3: red - 3
# B4: near infrared (nir) - 4
dev.off () # close others plots
# plotRGB function RGB colors
plotRGB (p224r63_2011, r=3, g=2, b=1, stretch="Lin")
# stretching as much possible the colors
# substitute red component with our NIR # in this way you can discriminete better the vegetation in the images, dark green is the forest, red and pink in general the pitch
plotRGB (p224r63_2011, r=4, g=3, b=2, stretch="Lin")
# EX. NIR on top of the G component # NIR false colours
plotRGB (p224r63_2011, r=3, g=4, b=2, stretch="Lin")
# dark forest is forest with high amount of water
# yellow area are bare in this way to set the colors
plotRGB (p224r63_2011, r=3, g=2, b=4, stretch="Lin")
################ day 2
#setwd("~/lab/") # linux
setwd("C:/lab/") # windows
# setwd("/Users/nome/Desktop/lab") # mac
load("R_code_remote_sensing.RData")
library(raster) # list
library (RStoolbox)
ls()
p224r63_1988 <- brick("p224r63_1988_masked.grd")
plot(p224r63_1988)
p224r63_2011 <- brick("p224r63_2011_masked.grd")
plot(p22r63_2011)
# Exercise: plot the image using the nir on the "r" component in the RGB space
plotRGB(p224r63_1988, r=4, g=3, b=2, stretch="Lin")
# Exrecise: plot in visible RGB 321 both images
par(mfrow=c(2,1))
plotRGB(p224r63_1988, r=3, g=2, b=1, stretch="Lin")
plotRGB(p224r63_2011, r=3, g=2, b=1, stretch="Lin")
# plot together two images, 1988 e 2011
par(mfrow=c(2,1))
plotRGB(p224r63_1988, r=4, g=3, b=2, stretch="Lin", main="1988")
plotRGB(p224r63_2011, r=4, g=3, b=2, stretch="Lin", main="2011")
# plot it in the way to see water, so with different colors relate at the NIR= water= 0
par(mfrow=c(2,1))
plotRGB(p224r63_1988, r=4, g=3, b=2, stretch="Lin")
plotRGB(p224r63_2011, r=4, g=3, b=2, stretch="Lin")
# enhance the noise!
par(mfrow=c(2,1))
plotRGB(p224r63_1988, r=4, g=3, b=2, stretch="hist")
plotRGB(p224r63_2011, r=4, g=3, b=2, stretch="hist")
# spectral indices
# dvi1988 = nir1988-red1988
dvi1988 <- p224r63_1988$B4_sre - p224r63_1988$B3_sre
plot(dvi1988)
# Exercise: calculate dvi for 2011
dvi2011 <- p224r63_2011$B4_sre - p224r63_2011$B3_sre
plot(dvi2011)
dev.off()
#giving colors at the div looking for the best grafic information
cl <- colorRampPalette(c("darkorchid3","light blue","lightpink4"))(100)
plot(dvi2011, col=cl)
#or
cldvi <- colorRampPalette(c("light blue","light green","green"))(100) # not the best for my purpose
plot(dvi2011, col=cldvi)
#or
cl2 <- colorRampPalette(c("yellow",'light blue','lightpink4'))(100)
plot(dvi2011,col=cl2)
# Exercise: dvi for 1988
par(mfrow= c(2,1))
plot(dvi1988, col=cl, ylab="1988")
plot(dvi2011, col=cl, ylab="2011")
par(mfrow= c(2,1))
# differances between 1988-2011
plot(dvi1988,col=cl3)
plot(dvi2011, col=cl)
dev.off()
diff <- dvi2011 - dvi1988
plot (diff)
# aggregate pixels
# changing the grain and resampling with factor=10 and 100 to minimize the weight
p224r63_2011res10 <- aggregate(p224r63_2011, fact=10)
p224r63_2011res100 <- aggregate(p224r63_2011, fact=100)
#EX. plot all together images with diff grain
par(mfrow=c(3,1))
plotRGB(p224r63_2011, r=4, g=3, b=2, stretch="Lin")
plotRGB(p224r63_2011res10, r=4, g=3, b=2, stretch="Lin")
plotRGB(p224r63_2011res100, r=4, g=3, b=2, stretch="Lin")
#information about image p224r63_2011
p224r63_2011