Vegetation Indices#
MODIS NDVI#
Data download from Google Earth Engine#
The MOD13Q1 V6 product provides a Vegetation Index (VI) value at a per pixel basis. There are two primary vegetation layers. The first is the Normalized Difference Vegetation Index (NDVI) which is referred to as the continuity index to the existing National Oceanic and Atmospheric Administration-Advanced Very High Resolution Radiometer (NOAA-AVHRR) derived NDVI. The second vegetation layer is the Enhanced Vegetation Index (EVI) that minimizes canopy background variations and maintains sensitivity over dense vegetation conditions. The MODIS NDVI and EVI products are computed from atmospherically corrected bi-directional surface reflectances that have been masked for water, clouds, heavy aerosols, and cloud shadows.
To compute and download the mean anual MODIS NDVI data from google earth engine. Open the Google Earth Engine Code and paste the lines of code provided below
1 /**** Start of imports. If edited, may not auto-convert in the playground. ****/
2 var table = ee.FeatureCollection("projects/ee-briansimiyu/assets/Africa_shapefile");
3 /***** End of imports. If edited, may not auto-convert in the playground. *****/
4 var dataset = ee.ImageCollection('MODIS/006/MOD13Q1')
5 .filter(ee.Filter.date('2018-01-01', '2019-01-01'))
6 .mean()
7 .clip(table);
8
9 var ndvi = dataset.select('NDVI');
10 var ndviVis = {
11 min: 0.0,
12 max: 8000.0,
13 palette: [
14 'FFFFFF', 'CE7E45', 'DF923D', 'F1B555', 'FCD163', '99B718', '74A901',
15 '66A000', '529400', '3E8601', '207401', '056201', '004C00', '023B01',
16 '012E01', '011D01', '011301'
17 ],
18 };
19 Map.centerObject(table);
20 Map.addLayer(ndvi, ndviVis, 'NDVI');
21
22 Export.image.toCloudStorage({
23 image:ndvi,
24 description: 'NDVI',
25 maxPixels:1e13,
26 scale:250,
27 bucket:'oss_ldms_vi',
28 region:table
29 })
Landsat derived NDVI, MSAVI, SAVI#
Data download from Google Earth Engine#
Several annual NDVI composite techniques have been discussed to overcome current Landsat 5 artifacts. MISLAND uses annual NDVI products based on different percentiles in order to better qualify and quantify the Vegetation Loss Index.
To compute and download desired percentile composites from google earth engine. Open the Google Earth Engine Code and paste the lines of code provided below
1 /// Required Data Inputs
2 // ===================
3 // * USGS/NASA's Landsat 4 surface reflectance tier 1 dataset (August 1982 - December 1993)
4 // * USGS/NASA's Landsat 5 surface reflectance tier 1 dataset (January 1, 1984 - May 5, 2012)
5 // * USGS/NASA's Landsat 7 surface reflectance tier 1 dataset (January 1, 1999 - December 31, 2019)
6 // * USGS/NASA's Landsat 8 surface reflectance tier 1 dataset (April 11, 2013 - December 31, 2019)
7 // * Study Area Polygon
8
9 var countries = ee.FeatureCollection("USDOS/LSIB_SIMPLE/2017");
10
11 var Year='2000'
12
13
14 var country = 'Libya';
15 var ALGO = "PERCENTILE65"; // PERCENTILE75 // PERCENTILE65 // PERCENTILE60 // MEDIAN
16 var cloudCoveragePercentage = 80;
17
18
19 var studyArea = countries.filter(ee.Filter.eq('country_na',country ))
20 //Map.addLayer(studyArea);
21
22 /*
23 borders are quite coarse
24 var northAfrica = ee.FeatureCollection('users/derickongeri/Admin')
25 var country = 'TUNISIA';
26 var studyArea = northAfrica.filter(ee.Filter.eq('NAME', country))
27 Map.addLayer(studyArea);
28 */
29
30
31
32
33 var start_date = Year+ '-01-01';
34 var end_date = Year+ '-12-31';
35
36
37 //--------------------------------------------------------------------
38 // Landsat 4, 5, 7 cloudmask
39 //--------------------------------------------------------------------
40
41 // If the cloud bit (5) is set and the cloud confidence (7) is high
42 // or the cloud shadow bit is set (3), then it's a bad pixel.
43 var cloudMaskL7 = function(image) {
44 var qa = image.select('QA_PIXEL');
45 var cloud = qa.bitwiseAnd(1 << 5)
46 .and(qa.bitwiseAnd(1 << 7))
47 .or(qa.bitwiseAnd(1 << 3));
48
49 // Remove edge pixels that don't occur in all bands
50 //var mask2 = image.mask().reduce(ee.Reducer.min())//.focal_min(300,'square','meters').eq(0);
51 //var mask2 = image.select('B4').reduce(ee.Reducer.min()).gt(0)//.focal_min(500,'square','meters');
52 // Remove edge pixels that don't occur in all bands
53 var mask3 =
54 (image.select('B3').gt(100))
55 .and(image.select('B4').gt(100))
56
57
58 .and(image.select('B4').lt(10000))
59 .and(image.select('B3').lt(10000))
60
61
62
63 return image.updateMask(cloud.not()).updateMask(mask3)//.updateMask(mask2)//.clip(image.geometry().buffer(-5000))//.or(mask3));
64 };
65
66
67 var cloudMaskL45 = function(image) {
68 var qa = image.select('QA_PIXEL');
69 var cloud = qa.bitwiseAnd(1 << 5)
70 .and(qa.bitwiseAnd(1 << 7))
71 .or(qa.bitwiseAnd(1 << 3));
72
73 // Remove edge pixels that don't occur in all bands
74 //var mask2 = image.mask().reduce(ee.Reducer.min());
75 var mask2 =
76 (image.select('B3').gt(100))
77 .and(image.select('B4').gt(100))
78
79
80 .and(image.select('B4').lt(10000))
81 .and(image.select('B3').lt(10000))
82
83 return (image.updateMask(cloud.not()).updateMask(mask2))//.clip(image.geometry().buffer(-5000))//.updateMask(mask2);
84 };
85
86
87 //--------------------------------------------------------------------
88 // Landsat 8 cloudmask
89 //--------------------------------------------------------------------
90
91 // Bits 3 and 5 are cloud shadow and cloud, respectively.
92 function maskL8sr(image) {
93 var cloudShadowBitMask = (1 << 3);
94 var cloudsBitMask = (1 << 5);
95
96 // Get the pixel QA band.
97 var qa = image.select('pixel_qa');
98
99 // Both flags should be set to zero, indicating clear conditions.
100 var mask = qa.bitwiseAnd(cloudShadowBitMask).eq(0)
101 .and(qa.bitwiseAnd(cloudsBitMask).eq(0));
102 var mask2 =
103
104 (image.select('B5').gt(100))
105 .and(image.select('B4').gt(100))
106
107
108 .and(image.select('B5').lt(10000))
109 .and(image.select('B4').lt(10000))
110
111 //var mask2 = image.mask().reduce(ee.Reducer.min()).focal_min(500,'square','meters');
112 //return image
113 return image.updateMask(mask).updateMask(mask2)//.clip(image.geometry().buffer(-5000));
114 }
115
116
117
118
119 // Apply Cloudmask to L4.5.7
120 var L4 = ee.ImageCollection("LANDSAT/LT04/C02/T2_TOA")
121 .filterDate(start_date, end_date)
122 .filter(ee.Filter.lessThan('CLOUD_COVER_LAND', cloudCoveragePercentage))
123 .filterBounds(studyArea)
124 .map(cloudMaskL45)
125 .select(['B3', 'B4'], ['RED', 'NIR']);;
126
127 var L5 = ee.ImageCollection('LANDSAT/LT05/C02/T2_TOA')
128 .filterDate(start_date, end_date)
129 .filter(ee.Filter.lessThan('CLOUD_COVER_LAND', cloudCoveragePercentage))
130 .filterBounds(studyArea)
131 .map(cloudMaskL45)
132 .select(['B3', 'B4'], ['RED', 'NIR']);;
133
134 var L7a = ee.ImageCollection('LANDSAT/LE07/C02/T2_TOA')
135 .filterDate('1999-01-01', '2003-04-01')
136 .filterDate(start_date, end_date)
137 .filter(ee.Filter.lessThan('CLOUD_COVER_LAND', 100))
138 .filterBounds(studyArea)
139 .map(cloudMaskL7)
140 .select(['B3', 'B4'], ['RED', 'NIR']);;
141 var L7b = ee.ImageCollection('LANDSAT/LE07/C02/T2_TOA')
142 .filterDate('2012-01-01', '2013-12-31')
143 .filterDate(start_date, end_date)
144 .filter(ee.Filter.lessThan('CLOUD_COVER_LAND', 100))
145 .filterBounds(studyArea)
146 .map(cloudMaskL7)
147 .select(['B3', 'B4'], ['RED', 'NIR']);;
148
149 var L7 = L7a.merge(L7b);
150
151 var L8 = ee.ImageCollection('LANDSAT/LC08/C02/T2_TOA')
152 .filterDate(start_date, end_date)
153 .filter(ee.Filter.lessThan('CLOUD_COVER', cloudCoveragePercentage))
154 .filterBounds(studyArea)
155 //.filterBounds(AOI)
156 .map(maskL8sr)
157 .select(['B4', 'B5'], ['RED', 'NIR']);
158
159
160
161
162
163
164 //Define collection
165
166
167 //--------------------------------------------------------------------
168 // Merge Landsat 4, 5, 8 imagery collections and filter all by date/place
169 //--------------------------------------------------------------------
170
171 //Merge Landsat 4, 5 , 7 '
172
173 var L4578 = L4.merge(L5).merge(L7).merge(L8);
174 print(L4578)
175
176
177 //--------------------------------------------------------------------
178 // Create NDVI Collection
179 //--------------------------------------------------------------------
180
181 var NDVI = function(image) {
182 return image.normalizedDifference(['NIR', 'RED']).rename('NDVI');
183 //return image.addBands(ndvi);
184 };
185
186 if (ALGO=='MEDIAN'){
187 var suffix = 'median';
188 var annualNDVI = L4578.map(NDVI).median().clip(studyArea);
189 }
190
191 if (ALGO=='PERCENTILE75'){
192 var suffix = '75pc';
193 var annualNDVI = L4578.map(NDVI).reduce(ee.Reducer.percentile([75])).clip(studyArea);
194 }
195
196 if (ALGO=='PERCENTILE65'){
197 var suffix = '65pc';
198 var annualNDVI = L4578.map(NDVI).reduce(ee.Reducer.percentile([65])).clip(studyArea);
199 }
200
201
202 var ndvi_visualization = {
203 min: -0.22789797020331423,
204 max: 0.6575894075894075,
205 palette: 'FFFFFF, CE7E45, DF923D, F1B555, FCD163, 99B718, 74A901, 66A000, 529400,' +
206 '3E8601, 207401, 056201, 004C00, 023B01, 012E01, 011D01, 011301'
207 };
208 Map.addLayer(annualNDVI, ndvi_visualization, 'NDVI');
209
210 //--------------------------------------------------------------------
211 // Export as GeoTIFF
212 //--------------------------------------------------------------------
213
214
215 Export.image.toDrive({
216 image: annualNDVI,
217 description: country + '_NDVI_' + suffix + '_' + Year,
218 scale: 30,
219 region: studyArea,
220 maxPixels: 1e13,
221 fileFormat: 'GeoTIFF',
222 folder:'GEE_classification',
223 formatOptions: {
224 cloudOptimized: true
225 },
226 skipEmptyTiles: true
227 });
228
229 //Map.addLayer(L7.first(), {}, 'L7');