-
Notifications
You must be signed in to change notification settings - Fork 0
/
transformations
428 lines (388 loc) · 13.1 KB
/
transformations
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
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
345
346
347
348
349
350
351
352
353
354
355
356
357
358
359
360
361
362
363
364
365
366
367
368
369
370
371
372
373
374
375
376
377
378
379
380
381
382
383
384
385
386
387
388
389
390
391
392
393
394
395
396
397
398
399
400
401
402
403
404
405
406
407
408
409
410
411
412
413
414
415
416
417
418
419
420
421
422
423
424
425
426
427
428
//Advanced Pixel-Based Image
//Transformations
//////////////////////The Enhanced Vegetation Index (EVI) calculation
/*to examine variation
in the vegetation in these regions*/
// Import and filter imagery by location and date.
var sfoPoint = ee.Geometry.Point(-122.3774, 37.6194);
var sfoImage = ee.ImageCollection('COPERNICUS/S2')
.filterBounds(sfoPoint)
.filterDate('2020-02-01', '2020-04-01')
.first();
Map.centerObject(sfoImage, 11);
// Calculate EVI using Sentinel 2
// Extract the bands and divide by 10,000 to account for scaling done.
var nirScaled = sfoImage.select('B8').divide(10000);
var redScaled = sfoImage.select('B4').divide(10000);
var blueScaled = sfoImage.select('B2').divide(10000);
// Calculate the numerator, note that order goes from left to right.
var numeratorEVI =
(nirScaled.subtract(redScaled)).multiply(2.5);
// Calculate the denominator.
var denomClause1 = redScaled.multiply(6);
var denomClause2 = blueScaled.multiply(7.5);
var denominatorEVI = nirScaled.add(denomClause1)
.subtract(denomClause2).add(1);
// Calculate EVI and name it.
var EVI =
numeratorEVI.divide(denominatorEVI).rename('EVI');
// And now map EVI using our vegetation palette.
var vegPalette = ['red', 'white', 'green'];
var visParams = {min: -1, max: 1, palette: vegPalette};
Map.addLayer(EVI, visParams, 'EVI');
// Calculate EVI.
var eviExpression = sfoImage.expression({
expression: '2.5 * ((NIR -RED) / (NIR + 6 * RED -7.5 * BLUE + 1))',
map: { // Map between variables in the expression and images.
'NIR': sfoImage.select('B8').divide(10000),
'RED': sfoImage.select('B4').divide(10000),
'BLUE': sfoImage.select('B2').divide(10000)
}
});
// And now map EVI using our vegetation palette.
Map.addLayer(eviExpression, visParams, 'EVI Expression');
/////////////////////Burned Area Index (BAI)
/*
To examine burn indices, load an image from 2013 showing the Rim Fire in the
Sierra Nevada, California mountains. We will use Landsat 8 to explore this fire.
*/
// Examine the true-color Landsat 8 images for the 2013 Rim Fire.
var burnImage =
ee.ImageCollection('LANDSAT/LC08/C02/T1_TOA')
.filterBounds(ee.Geometry.Point(-120.083, 37.850))
.filterDate('2013-09-15', '2013-09-27')
.sort('CLOUD_COVER')
.first();
Map.centerObject(ee.Geometry.Point(-120.083, 37.850), 11);
var rgbParams = {
bands: ['B4', 'B3', 'B2'],
min: 0,
max: 0.3
};
Map.addLayer(burnImage, rgbParams, 'True-Color Burn Image');
// Calculate BAI.
var bai = burnImage.expression(
'1.0 / ((0.1 -RED)**2 + (0.06 -NIR)**2)', {
'NIR': burnImage.select('B5'),
'RED': burnImage.select('B4'),
});
// Display the BAI image.
var burnPalette = ['green', 'blue', 'yellow', 'red'];
Map.addLayer(bai, {
min: 0,
max: 400,
palette: burnPalette
}, 'BAI');
// Manipulating Images with Matrix Algebra
/////
// Manipulating images with matrices
/////
// Begin Tasseled Cap example.
var landsat5RT = ee.Array([
[0.3037, 0.2793, 0.4743, 0.5585, 0.5082, 0.1863],
[-0.2848, -0.2435, -0.5436, 0.7243, 0.0840, -0.1800],
[0.1509, 0.1973, 0.3279, 0.3406, -0.7112, -0.4572],
[-0.8242, 0.0849, 0.4392, -0.0580, 0.2012, -0.2768],
[-0.3280, 0.0549, 0.1075, 0.1855, -0.4357, 0.8085],
[0.1084, -0.9022, 0.4120, 0.0573, -0.0251, 0.0238]
]);
print('RT for Landsat 5', landsat5RT);
/************ tasseled cap transformation*/
//in Odessa, WA, USA
// Define a point of interest in Odessa, Washington, USA.
var point = ee.Geometry.Point([-118.7436019417829,
47.18135755009023]);
Map.centerObject(point, 10);
// Filter to get a cloud free image to use for the TC.
var imageL5 = ee.ImageCollection('LANDSAT/LT05/C02/T1_TOA')
.filterBounds(point)
.filterDate('2008-06-01', '2008-09-01')
.sort('CLOUD_COVER')
.first();
//Display the true-color image.
var trueColor = {
bands: ['B3', 'B2', 'B1'],
min: 0,
max: 0.3
};
Map.addLayer(imageL5, trueColor, 'L5 true color');
var bands = ['B1', 'B2', 'B3', 'B4', 'B5', 'B7'];
// Make an Array Image, with a one dimensional array per pixel.
// This is essentially a list of values of length 6,
// one from each band in variable 'bands.'
var arrayImage1D = imageL5.select(bands).toArray();
// Make an Array Image with a two dimensional array per
//pixel,
// of dimensions 6x1. This is essentially a one column matrix with
// six rows, with one value from each band in 'bands.'
// This step is needed for matrix multiplication (p0).
var arrayImage2D = arrayImage1D.toArray(1);
//Multiply RT by p0.
var tasselCapImage = ee.Image(landsat5RT)
// Multiply the tasseled cap coefficients by the array
// made from the 6 bands for each pixel.
.matrixMultiply(arrayImage2D)
// Get rid of the extra dimensions.
.arrayProject([0])
// Get a multi-band image with TC-named bands.
.arrayFlatten(
[
['brightness', 'greenness', 'wetness',
'fourth', 'fifth',
'sixth'
]
]);
var vizParams = {
bands: ['brightness', 'greenness', 'wetness'],
min: -0.1,
max: [0.5, 0.1, 0.1]
};
Map.addLayer(tasselCapImage, vizParams, 'TC components');
///////////////////////////////////
//// principal component analysis (PCA)
//////////////////////////////////
// Begin PCA example.
// Select and map a true-color L8 image.
var imageL8 = ee.ImageCollection('LANDSAT/LC08/C02/T1_TOA')
.filterBounds(point)
.filterDate('2018-06-01', '2018-09-01')
.sort('CLOUD_COVER')
.first();
var trueColorL8 = {
bands: ['B4', 'B3', 'B2'],
min: 0,
max: 0.3
};
Map.addLayer(imageL8, trueColorL8, 'L8 true color');
// Select which bands to use for the PCA.
var PCAbands = ['B2', 'B3', 'B4', 'B5', 'B6', 'B7', 'B10',
'B11'];
// Convert the Landsat 8 image to a 2D array for the later matrix
// computations.
var arrayImage = imageL8.select(PCAbands).toArray();
// Calculate the covariance using the reduceRegion method.
var covar = arrayImage.reduceRegion({
reducer: ee.Reducer.covariance(),
maxPixels: 1e9
});
// Extract the covariance matrix and store it as an array.
var covarArray = ee.Array(covar.get('array'));
//Compute and extract the eigenvectors
var eigens = covarArray.eigen();
var eigenVectors = eigens.slice(1, 1);
// Perform matrix multiplication
var principalComponents = ee.Image(eigenVectors)
.matrixMultiply(arrayImage.toArray(1));
var pcImage = principalComponents
// Throw out an unneeded dimension, [[]] -> [].
.arrayProject([0])
// Make the one band array image a multi-band image, [] -> image.
.arrayFlatten([
['pc1', 'pc2', 'pc3', 'pc4', 'pc5', 'pc6', 'pc7',
'pc8']
]);
// Stretch this to the appropriate scale.
Map.addLayer(pcImage.select('pc1'), {}, 'pc1');
var visParamsPCA = {
bands: ['pc1', 'pc3', 'pc4'],
min: [-455.09, -2.206, -4.53],
max: [-417.59, -1.3, -4.18]
};
Map.addLayer(pcImage, visParamsPCA, 'PC_multi');
//////////////////////////////////////
// Create and print a uniform kernel to see its weights.
print('A uniform kernel:', ee.Kernel.square(2));
// Define a point of interest in Odessa, Washington, USA.
var point = ee.Geometry.Point([-118.71845096212049,
47.15743083101999]);
Map.centerObject(point);
// Load NAIP data.
var imageNAIP = ee.ImageCollection('USDA/NAIP/DOQQ')
.filterBounds(point)
.filter(ee.Filter.date('2017-01-01', '2018-12-31'))
.first();
Map.centerObject(point, 17);
var trueColor = {
bands: ['R', 'G', 'B'],
min: 0,
max: 255
};
Map.addLayer(imageNAIP, trueColor, 'true color');
// Begin smoothing example.
// Define a square, uniform kernel.
var uniformKernel = ee.Kernel.square({
radius: 2,
units: 'meters', });
// Convolve the image by convolving with the smoothing kernel.
var smoothed = imageNAIP.convolve(uniformKernel);
Map.addLayer(smoothed, {
min: 0,
max: 255
}, 'smoothed image');
///////////////////////Gaussian Smoothing
// Begin Gaussian smoothing example.
// Print a Gaussian kernel to see its weights.
print('A Gaussian kernel:', ee.Kernel.gaussian(2));
// Define a square Gaussian kernel:
var gaussianKernel = ee.Kernel.gaussian({
radius: 2,
units: 'meters',
});
// Convolve the image with the Gaussian kernel.
var gaussian = imageNAIP.convolve(gaussianKernel);
Map.addLayer(gaussian, {
min: 0,
max: 255
}, 'Gaussian smoothed image');
/////////// Edge Detection
// Begin edge detection example.
// For edge detection, define a Laplacian kernel.
var laplacianKernel = ee.Kernel.laplacian8();
// Print the kernel to see its weights.
print('Edge detection Laplacian kernel:', laplacianKernel);
// Convolve the image with the Laplacian kernel.
var edges = imageNAIP.convolve(laplacianKernel);
Map.addLayer(edges, {
min: 0,
max: 255
}, 'Laplacian convolution image');
// Begin image sharpening example.
// Define a "fat" Gaussian kernel.
var fat = ee.Kernel.gaussian({
radius: 3,
sigma: 3,
magnitude: -1,
units: 'meters'
});
// Define a "skinny" Gaussian kernel.
var skinny = ee.Kernel.gaussian({
radius: 3,
sigma: 0.5,
units: 'meters' });
// Compute a difference-of-Gaussians (DOG) kernel.
var dog = fat.add(skinny);
// Print the kernel to see its weights.
print('DoG kernel for image sharpening',dog);
// Add the DoG convolved image to the original image.
var sharpened = imageNAIP.add(imageNAIP.convolve(dog));
Map.addLayer(sharpened, {
min: 0,
max: 255
}, 'DoG edge enhancement');
/////////////////////////////////
////////////Nonlinear Convolution
// Begin median example.
// Pass a median neighborhood filter using our uniformKernel.
var median = imageNAIP.reduceNeighborhood({
reducer: ee.Reducer.median(),
kernel: uniformKernel
});
Map.addLayer(median, {
min: 0,
max: 255
}, 'Median Neighborhood Filter');
// Mode example
// Create and display a simple two-class image.
var veg = imageNAIP.select('N').gt(200);
// Display the two-class (binary) result.
var binaryVis = {
min: 0,
max: 1,
palette: ['black', 'green']
};
Map.addLayer(veg, binaryVis, 'Vegetation categorical image');
// Compute the mode in each 5x5 neighborhood and display the result.
var mode = veg.reduceNeighborhood({
reducer: ee.Reducer.mode(),
kernel: uniformKernel
});
Map.addLayer(mode, binaryVis, 'Mode Neighborhood Filter on Vegetation categorical image');
/////////////////Morphological Processing
// Begin Dilation example.
// Dilate by taking the max in each 5x5 neighborhood.
var max = veg.reduceNeighborhood({
reducer: ee.Reducer.max(),
kernel: uniformKernel
});
Map.addLayer(max, binaryVis, 'Dilation using max');
// Begin Erosion example.
// Erode by taking the min in each 5x5 neighborhood.
var min = veg.reduceNeighborhood({
reducer: ee.Reducer.min(),
kernel: uniformKernel
});
Map.addLayer(min, binaryVis, 'Erosion using min');
// Begin Opening example.
// Perform an opening by dilating the eroded image.
var openedVeg = min.reduceNeighborhood({
reducer: ee.Reducer.max(),
kernel: uniformKernel
});
Map.addLayer(openedVeg, binaryVis, 'Opened image');
// Begin Closing example.
// Perform a closing by eroding the dilated image.
var closedVeg = max.reduceNeighborhood({
reducer: ee.Reducer.min(),
kernel: uniformKernel
});
Map.addLayer(closedVeg, binaryVis, 'Closed image');
/////////////////Standard Deviation
// Begin Standard Deviation example.
// Define a big neighborhood with a 7-meter radius kernel.
var bigKernel = ee.Kernel.square({
radius: 7,
units: 'meters'
});
// Compute SD in a neighborhood.
var sd = imageNAIP.reduceNeighborhood({
reducer: ee.Reducer.stdDev(),
kernel: bigKernel
});
Map.addLayer(sd, {
min: 0,
max: 70
}, 'SD');
// Begin entropy example.
// Create an integer version of the NAIP image.
var intNAIP = imageNAIP.int();
// Compute entropy in a neighborhood.
var entropy = intNAIP.select('N').entropy(bigKernel);
Map.addLayer(entropy, {
min: 1,
max: 3
}, 'entropy');
// Begin GLCM example.
// Use the GLCM to compute a large number of texture measures.
var glcmTexture = intNAIP.glcmTexture(7);
print('view the glcmTexture output',glcmTexture);
var contrastVis = {
bands: ['R_contrast', 'G_contrast', 'B_contrast'],
min: 40,
max: 1000
};
Map.addLayer(glcmTexture, contrastVis, 'contrast');
// Begin spatial statistics example using Geary's C.
// Create a list of weights for a 9x9 kernel.
var list = [1, 1, 1, 1, 1, 1, 1, 1, 1];
// The center of the kernel is zero.
var centerList = [1, 1, 1, 1, 0, 1, 1, 1, 1];
// Assemble a list of lists: the 9x9 kernel weights as a 2-D matrix.
var lists = [list, list, list, list, centerList, list,
list, list,list
];
// Create the kernel from the weights.
// Non-zero weights represent the spatial neighborhood.
var kernel = ee.Kernel.fixed(9, 9, lists, -4, -4, false);
// Use the max among bands as the input.
var maxBands = imageNAIP.reduce(ee.Reducer.max());
// Convert the neighborhood into multiple bands.
var neighBands = maxBands.neighborhoodToBands(kernel);
// Compute local Geary's C, a measure of spatial association.
var gearys =
maxBands.subtract(neighBands).pow(2).reduce(ee.Reducer
.sum())
.divide(Math.pow(9, 2));
Map.addLayer(gearys, {
min: 20,
max: 2500
}, "Geary's C");