Analysis of change detection of Jharia CoalFields using Google Earth Engine
Jharia is a notified area and one of eight development blocks of Dhanbad district in Jharkhand state, India. Jharia is famous for its rich coal resources, used to make coke. The coalfield lies in the Damodar River Valley covering about 110 square miles (280 square km) and produces bituminous coal suitable for coke and most of India’s coal comes from Jharia. But its surface and subsurface coal fire is the major problem here like in many coal-producing countries such as China, Australia, India, etc.
In India, at Jharia some of the subsurface coal fires have already been burning over many decades, and now they are spreading caused by mining introduced spontaneous combustion. Coal fire is the main cause of Jharia that has burned underground coal for nearly a century. The first fire was detected in
1916. In 1972, more than 70 mine fires were reported in this region and as reported in 2010, 68 fires were burning beneath a 58-square mile (150 km²) region of the Jharia coalfield (JCF), this field was the only major source of coke coal. Geographical information systems (GIS) and remote sensing are well-established Information Technologies, whose applications in land and natural resource management are widely recognized. These technologies provide a cost-effective and accurate alternative to conventional approaches used for understanding landscape dynamics.
Digital change detection techniques based on multi-temporal and multi-spectral remotely sensed data have demonstrated great potential as a means to understanding landscape dynamics to detect, identify, map, and monitor differences in land use/land cover patterns over time, irrespective of the causal factors. Recent improvements in satellite image quality and availability have made it possible to perform image analysis at a much larger scale than in the past. A detailed understanding of the impact of mining on changes in land use/land cover pattern has become necessary for the Jharia Coal Field. The present study was undertaken to analyze the extent of human-induced landscape transformation in the mining-affected areas of Jharia Coal Field by interpreting temporal remote sensing data.
Region of Interest (study Area)
Jharia Coal Field area has attracted people over the past centuries due to its large amount of bituminous coal(suitable for coke) production in India.
The sickle-shaped Jharia coalfield (JCF) is located in Jharkhand’s Dhanbad district, roughly 250 kilometers northwest of Kolkata. with an area of 700 km², the JCF has bordered by 23 0 37'N to 23 0 50’N latitude and 86⁰ 8' E to 86⁰ 30' E longitudes.
The major part of the coalfield (about 400 km²) lies in the Dhanbad district of Jharkhand. Coalfield is connected by Major Highways Road with Ranchi (117 km), Asansol (60 km), Jamshedpur (108 km), and Dhanbad (8 km).
Google Earth Engine
Google Earth Engine is a computing platform that allows users to run geospatial analysis on Google’s infrastructure. There are several ways to interact with the platform. The Code Editor is a web-based IDE for writing and running scripts. The Explorer is a lightweight web app for exploring our data catalog and running simple analyses. The client libraries provide Python and JavaScript wrappers around our web API.Google Earth Engine is a geospatial processing service. With Earth Engine, you can perform geospatial processing at scale, powered by the Google Cloud Platform. The purpose of Earth Engine is to:
- Provide an interactive platform for geospatial algorithm development at scale
- Enable high-impact, data-driven science
- Make substantive progress on global challenges that involve a large geospatial dataset
LANDSAT 7
Landsat 7 was launched from Vandenberg Air Force Base in California on April 15, 1999, on a Delta II rocket. The satellite carries the Enhanced Thematic Mapper Plus (ETM+) sensor. This instrument was improved from previous instrumentation designs. The primary features on Landsat 7 include a panchromatic band with 15-meter spatial resolution, an onboard full aperture solar calibrator, five percent absolute radiometric calibration, and a
thermal infrared channel with a four-fold improvement in spatial resolution over Thematic Mapper (TM).
Landsat 7 Satellite Orbit Facts
- Orbits the Earth at 705 km (438 mi) in a sun-synchronous, near-polar orbit (98.2⁰ inclination)
- Circles the Earth every 99 minutes
- Has a 16-day repeat cycle with an equatorial crossing time: 10:00 a.m. (+/- 15 minutes)
- Acquired on the Worldwide Reference System-2 (WRS-2) path/row system, with swath overlap (or sidelap) varying
from 7 percent at the Equator to a maximum of approximately 85 percent at extreme latitudes
NDVI
Normalized Difference Vegetation Index (NDVI) quantifies vegetation by measuring the difference between near-infrared (which vegetation strongly reflects) and red light (which vegetation absorbs). NDVI always ranges from -1 to+1. But there isn’t a distinct boundary for each type of land cover. NDVI is the most common index that analysts use in remote sensing.
For example, when you have negative values, it’s highly likely that it’s water. On the other hand, if you have an NDVI value close to +1, there’s a high possibility that it’s dense green leaves. But when NDVI is close to zero, there aren’t green leaves and it could even be an urbanized area. Normalized Difference Vegetation Index (NDVI) uses the NIR and red channels in its formula.
NDVI = (NIR — Red) / (NIR + Red)Healthy vegetation (chlorophyll) reflects more near-infrared (NIR) and green light compared to other wavelengths. But it absorbs more red and blue light.
This is why our eyes see vegetation as green color. If you could see near-infrared, then it would be strong for vegetation too. Satellite sensors
like Landsat and Sentinel-2 both have the necessary bands with NIR and red.
The result of this formula generates a value between -1 and +1. If you have low reflectance (or low values) in the red channel and high reflectance in the NIR channel, this will yield a high NDVI value. And vice versa. Overall, NDVI is a standardized way to measure healthy vegetation. When you have high NDVI values, you have healthier vegetation. When you have low NDVI, you have less or no vegetation. Generally, if you want to see the vegetation change over time, then you will have to perform the atmospheric correction.
BU and NDBI
Build-up Index is the index for analysis of urban patterns using NDBI and NDVI. The built-up index is the binary image
with only a higher positive value indicating built-up and barren thus, allows BU to map the built-up area automatically.
- The build-up areas and bare soil reflect more SWIR than NIR.
- Waterbody doesn’t reflect on the Infrared spectrum.
- In the case of the greenie surface, a reflection of NIR is higher than the SWIR spectrum
BU = NDBI — NDVIBU = NDBI — NDVI
NDBI = (SWIR — NIR) / (SWIR + NIR)
For Landsat 7 data, NDBI = (Band 5 — Band 4) / (Band 5 + Band 4)
For Landsat 8 data, NDBI = (Band 6 — Band 5) / (Band 6 + Band 5)
Also, the Normalize Difference Build-up Index value lies between -1 to +1. A negative value of NDBI represents water
bodies whereas a higher value represents build-up areas. NDBI value for vegetation is low.
Bu_2000: 113924*(900 * 10-⁶ ) = 102.5316 sq. km
Bu_2005 : 156970*(900 * 10-⁶) = 141.2730 sq. km
Bu_2010 : 144422*(900 * 10-⁶) = 129.9798 sq. km
Bu_2015 : 156906*(900 * 10-⁶) = 141.2154 sq. km
Bu_2019 : 134454*(900 * 10-⁶) = 121.0086 sq. km
Total area of study = 173.6996 sq. km
Conclusion
These land-use/ land-cover maps provide information about the land use pattern of the study area. Five land use classes have been identified as built-up, water bodies, Vegetation, Dense Vegetation & bare soil. Change analysis shows that the built-up area has increased by 11.51%, the vegetation cover has decreased by 5.50% and water bodies have decreased by 28.57%. Poor planning of overburden disposal of the large open-cast mines contributed adversely to the increase of barren land & overburden area.
Jharia Coal Field is the third-largest coking coal-producing area in the world. It is experiencing rapid urbanization as well as mining activities as an increase in a built-up area and barren land & overburden in the last two decades. Thus, Information on land-use/land cover and possibilities for their optimal use is essential for the selection, planning, and implementation of land use schemes to meet the increasing demands for basic human needs and welfare in the state of increasing mining activities in this area.
var geometry = /* color: #23d63d */ee.Geometry.Polygon(
[[[86.20961922282639, 23.80437862483123],
[86.18592995280686, 23.792127527682606],
[86.1718537198967, 23.785844465688715],
[86.17048042888108, 23.777047668699673],
[86.17322701091233, 23.772020660373883],
[86.17837685222092, 23.769192882811495],
[86.18524330729905, 23.766050835652667],
[86.19176643962327, 23.76510820670878],
[86.20275276774827, 23.767936073053605],
[86.2264420377678, 23.77516256334672],
[86.24635475749436, 23.779875275494074],
[86.2594010221428, 23.784901980194256],
[86.27382057780686, 23.786158626001704],
[86.28034371013108, 23.785530304616458],
[86.29133003825608, 23.78301698870662],
[86.30128639811936, 23.778932746727563],
[86.31124275798264, 23.77264904703996],
[86.3201691495842, 23.770135482162125],
[86.33458870524827, 23.771706465902483],
[86.34660500163498, 23.770763877934908],
[86.3586212980217, 23.770135482162125],
[86.37819069499436, 23.764165570936246],
[86.38265389079514, 23.757252699970227],
[86.3857437955803, 23.746254194085534],
[86.38814705485764, 23.73965464466331],
[86.39535683268967, 23.731483310744732],
[86.39810341472092, 23.725825933270254],
[86.3970734464592, 23.713567439854337],
[86.38780373210373, 23.703822687292067],
[86.3805939542717, 23.694311636032914],
[86.37990730876389, 23.684251040410352],
[86.38643044108811, 23.675447383539286],
[86.39741676921311, 23.67135976978636],
[86.4080597745842, 23.667900919833876],
[86.42350929850998, 23.659410627403325],
[86.43071907634202, 23.651863238016013],
[86.44857185954514, 23.64777488694244],
[86.4629914152092, 23.648089380022103],
[86.46505135173264, 23.654379082866978],
[86.45784157390061, 23.663498614802002],
[86.46230476970139, 23.66947313570403],
[86.45646828288498, 23.678277195084906],
[86.44857185954514, 23.686766261976388],
[86.44788521403733, 23.695569155985314],
[86.45612496013108, 23.705314524878194],
[86.45955818767014, 23.714430501457024],
[86.45372170085373, 23.735174831683537],
[86.44410866374436, 23.757487015414693],
[86.43483894938889, 23.770683972110536],
[86.42144936198655, 23.7841937552114],
[86.39947670573655, 23.797388002137293],
[86.37681740397873, 23.805241084481608],
[86.35965126628342, 23.80838218447629],
[86.33973854655686, 23.8115232084921],
[86.32497566813889, 23.815920514455286],
[86.31364601725998, 23.81560642611177],
[86.29956978434983, 23.813721880091695],
[86.2762238370842, 23.81215140417713],
[86.25322121257248, 23.80995270598203],
[86.23639839763108, 23.808068077895744],
[86.21991890544358, 23.805555197899903]]]);//Setting visualization parameters for Landsat-7 (l7) and
// Landsat-8 (l8) false color composite (fcc)
var l7_fcc_vis = {min: 0, max: 3000, bands: ['B4', 'B3', 'B2']};
var l8_fcc_vis = {min: 0, max: 3000, bands: ['B5', 'B4', 'B3']};//Setting up region of interest
// GODHAR LAT LOG 86.3957, 23.77954
var roi = ee.Geometry.Point([86.3957, 23.77954])
.buffer(5000);
Map.addLayer(roi, {}, "ROI", false);
Map.centerObject(roi, 12);
var empty = ee.Image().byte();
var outline = empty.paint({
featureCollection: geometry, //--- specify the user defined area variable here
color: 1,
width: 3
});//Filtering L7 data for year 2000
var L7 = ee.ImageCollection("LANDSAT/LE07/C01/T1_SR")
// .filterDate('2000-01-01', '2000-12-31')
.filterDate('2000-01-01', '2005-01-01')
.filterBounds(roi)
.filterMetadata('CLOUD_COVER', 'less_than', 30);//Generate median (cloud free) image
var L7_2000_median_image = ee.Image(L7.median());
Map.addLayer(L7_2000_median_image, l7_fcc_vis, 'L7_2000_median_image', false);//Filtering L7 data for year 2005
var L7_2005 = ee.ImageCollection("LANDSAT/LE07/C01/T1_SR")
// .filterDate('2005-01-01', '2005-12-31')
.filterDate('2005-01-01', '2010-01-01')
.filterBounds(roi)
.filterMetadata('CLOUD_COVER', 'less_than', 30);//Generate median (cloud free) image
var L7_2005_median_image = ee.Image(L7_2005.median());
Map.addLayer(L7_2005_median_image, l7_fcc_vis, 'L7_2005_median_image', false);//Filtering L7 data for year 2010
var L7_2010 = ee.ImageCollection("LANDSAT/LE07/C01/T1_SR")
// .filterDate('2010-01-01', '2010-12-31')
.filterDate('2010-01-01', '2015-01-01')
.filterBounds(roi)
.filterMetadata('CLOUD_COVER', 'less_than', 30);//Generate median (cloud free) image
var L7_2010_median_image = ee.Image(L7_2010.median());
Map.addLayer(L7_2010_median_image, l7_fcc_vis, 'L7_2010_median_image', false);//Filtering L7 data for year 2015
var L7_2015 = ee.ImageCollection("LANDSAT/LE07/C01/T1_SR")
// .filterDate('2015-01-01', '2015-12-31')
.filterDate('2015-01-01', '2019-01-01')
.filterBounds(roi)
.filterMetadata('CLOUD_COVER', 'less_than', 30);//Generate median (cloud free) image
var L7_2015_median_image = ee.Image(L7_2015.median());
Map.addLayer(L7_2015_median_image, l7_fcc_vis, 'L7_2015_median_image', false);//Filtering L8 data for year 2019
var L8 = ee.ImageCollection("LANDSAT/LC08/C01/T1_SR")
.filterDate('2019-01-01', '2019-12-31')
.filterBounds(roi)
.filterMetadata('CLOUD_COVER', 'less_than', 30);//Generate median (cloud free) image
var L8_2019_median_image = ee.Image(L8.median());
Map.addLayer(L8_2019_median_image, l8_fcc_vis, 'L8_2019_median_image', false);// Zha 2003 (IJRS) ...
//NDVI_C is (NIR - Red)/(NIR + Red)
//NDVI_B is 245 for all pixels NDVI_C > 0
//NDBI_C is (SWIR – NIR)/(SWIR + NIR)
//NDBI_B is 254 for all pixels NDBI_C > 0
//BU is NDBI_B - NDVI_B; will have 254 for BU and Barrenvar ndvi_c_2000 = L7_2000_median_image
.normalizedDifference(['B4', 'B3']);
var ndvi_b_2000 = ee.Image.constant(254)
.updateMask(ndvi_c_2000.gt(0));var ndbi_c_2000 = L7_2000_median_image
.normalizedDifference(['B5', 'B4']);
var ndbi_b_2000 = ee.Image.constant(254)
.updateMask(ndbi_c_2000.gt(0));var bu_2000 = ndbi_b_2000.subtract(ndvi_b_2000).rename("bu_2000");Map.addLayer(ndvi_c_2000, {}, 'ndvi_c_2000', false);
Map.addLayer(ndvi_b_2000, {}, 'ndvi_b_2000', false);
Map.addLayer(ndbi_c_2000, {}, 'ndbi_c_2000', false);
Map.addLayer(ndbi_b_2000, {}, 'ndbi_b_2000', false);Map.addLayer(bu_2000, {}, 'bu_2000');// for 2005
var ndvi_c_2005 = L7_2005_median_image
.normalizedDifference(['B4', 'B3']);
var ndvi_b_2005 = ee.Image.constant(254)
.updateMask(ndvi_c_2005.gt(0));var ndbi_c_2005 = L7_2005_median_image
.normalizedDifference(['B5', 'B4']);
var ndbi_b_2005 = ee.Image.constant(254)
.updateMask(ndbi_c_2005.gt(0));var bu_2005 = ndbi_b_2005.subtract(ndvi_b_2005).rename("bu_2005");Map.addLayer(ndvi_c_2005, {}, 'ndvi_c_2005', false);
Map.addLayer(ndvi_b_2005, {}, 'ndvi_b_2005', false);
Map.addLayer(ndbi_c_2005, {}, 'ndbi_c_2005', false);
Map.addLayer(ndbi_b_2005, {}, 'ndbi_b_2005', false);Map.addLayer(bu_2005, {}, 'bu_2005', false);// for 2010
var ndvi_c_2010 = L7_2010_median_image
.normalizedDifference(['B4', 'B3']);
var ndvi_b_2010 = ee.Image.constant(254)
.updateMask(ndvi_c_2010.gt(0));var ndbi_c_2010 = L7_2010_median_image
.normalizedDifference(['B5', 'B4']);
var ndbi_b_2010 = ee.Image.constant(254)
.updateMask(ndbi_c_2010.gt(0));var bu_2010 = ndbi_b_2010.subtract(ndvi_b_2010).rename("bu_2010");Map.addLayer(ndvi_c_2010, {}, 'ndvi_c_2010', false);
Map.addLayer(ndvi_b_2010, {}, 'ndvi_b_2010', false);
Map.addLayer(ndbi_c_2010, {}, 'ndbi_c_2010', false);
Map.addLayer(ndbi_b_2010, {}, 'ndbi_b_2010', false);Map.addLayer(bu_2010, {}, 'bu_2010', false);// for 2015
var ndvi_c_2015 = L7_2015_median_image
.normalizedDifference(['B4', 'B3']);
var ndvi_b_2015 = ee.Image.constant(254)
.updateMask(ndvi_c_2015.gt(0));var ndbi_c_2015 = L7_2015_median_image
.normalizedDifference(['B5', 'B4']);
var ndbi_b_2015 = ee.Image.constant(254)
.updateMask(ndbi_c_2015.gt(0));var bu_2015 = ndbi_b_2015.subtract(ndvi_b_2015).rename("bu_2015");Map.addLayer(ndvi_c_2015, {}, 'ndvi_c_2015', false);
Map.addLayer(ndvi_b_2015, {}, 'ndvi_b_2015', false);
Map.addLayer(ndbi_c_2015, {}, 'ndbi_c_2015', false);
Map.addLayer(ndbi_b_2015, {}, 'ndbi_b_2015', false);Map.addLayer(bu_2015, {}, 'bu_2015', false);///ly for 2019var ndvi_c_2019 = L8_2019_median_image
.normalizedDifference(['B5', 'B4']);
var ndvi_b_2019 = ee.Image.constant(254)
.updateMask(ndvi_c_2019.gt(0));var ndbi_c_2019 = L8_2019_median_image
.normalizedDifference(['B6', 'B5']);
var ndbi_b_2019 = ee.Image.constant(254)
.updateMask(ndbi_c_2019.gt(0));var bu_2019 = ndbi_b_2019.subtract(ndvi_b_2019).rename("bu_2019");Map.addLayer(ndvi_c_2019, {}, 'ndvi_c_2019', false);
Map.addLayer(ndvi_b_2019, {}, 'ndvi_b_2019', false);
Map.addLayer(ndbi_c_2019, {}, 'ndbi_c_2019', false);
Map.addLayer(ndbi_b_2019, {}, 'ndbi_b_2019', false);
Map.addLayer(bu_2019, {}, 'bu_2019');
Map.addLayer(outline, {palette: 'FF0000'}, 'Study Area',true);
print("GEOMETRY area in km^2: ",ee.Number(geometry.area()).divide(1e6));print(bu_2000.reduceRegion(ee.Reducer.count(), geometry, 30));
print(bu_2005.reduceRegion(ee.Reducer.count(), geometry, 30));
print(bu_2010.reduceRegion(ee.Reducer.count(), geometry, 30));
print(bu_2015.reduceRegion(ee.Reducer.count(), geometry, 30));
print(bu_2019.reduceRegion(ee.Reducer.count(), geometry, 30));// var reg_chart=ui.Chart.image.seriesByRegion(L7,
// geometry,
// ee.Reducer.mean(),
// 'B4',
// 30,
// 'SENSING_TIME',
// 'CLOUD_COVER')
// .setOptions({
// title:' analysing data sr',
// vAxis:{title: 'vertical title'}
// });
// print(reg_chart);//chart 2000
var chart_2000 = ui.Chart.image.series({
imageCollection: L7.select('B4','B3','B2','B5'),
region: geometry,
reducer: ee.Reducer.mean(),
scale: 30
})
.setOptions({
lineWidth: 1,
title: 'year 2000-2005',
interpolateNulls: true,
vAxis: {title: 'measurse of band values B2,B3,B4'},
hAxis: {title: 'Date', format: 'YYYY-MMM'}
})
print(chart_2000);//chart 2005
var chart_2005 = ui.Chart.image.series({
imageCollection: L7_2005.select('B4','B3','B2','B5'),
region: geometry,
reducer: ee.Reducer.mean(),
scale: 30
})
.setOptions({
lineWidth: 1,
title: 'year 2005-2010',
interpolateNulls: true,
vAxis: {title: 'measurse of band values B2,B3,B4'},
hAxis: {title: 'Date', format: 'YYYY-MMM'}
})
print(chart_2005);//chart 2010
var chart_2010 = ui.Chart.image.series({
imageCollection: L7_2010.select('B4','B3','B2','B5'),
region: geometry,
reducer: ee.Reducer.mean(),
scale: 30
})
.setOptions({
lineWidth: 1,
title: 'year 2010-2015',
interpolateNulls: true,
vAxis: {title: 'measurse of band values B2,B3,B4'},
hAxis: {title: 'Date', format: 'YYYY-MMM'}
})
print(chart_2010);//chart 2015
var chart_2015 = ui.Chart.image.series({
imageCollection: L7_2015.select('B4','B3','B2','B5'),
region: geometry,
reducer: ee.Reducer.mean(),
scale: 30
})
.setOptions({
lineWidth: 1,
title: 'year 2015-2019',
interpolateNulls: true,
vAxis: {title: 'measurse of band values B2,B3,B4'},
hAxis: {title: 'Date', format: 'YYYY-MMM'}
})
print(chart_2015);//chart 2019
var chart_2019 = ui.Chart.image.series({
imageCollection: L8.select('B3','B4','B5','B6'),
region: geometry,
reducer: ee.Reducer.mean(),
scale: 30
})
.setOptions({
lineWidth: 1,
title: 'year 2019',
interpolateNulls: true,
vAxis: {title: 'measurse of band values B2,B3,B4'},
hAxis: {title: 'Date', format: 'YYYY-MMM'}
})
print(chart_2019);//2000 = 126 km sq area
//2019 = 150 km sq area
Special thanks: Kunal, Ranjan, Rahul, Ranjeet, Sanju, Bimal &
Prof.Chinmay Mondal ,department of Mining Engineering BIT Sindri .
