Thursday, August 20, 2020

GEE: extarct and plot data from sentinel 2A_B imagery

 

First, Let us pull out a toa image from the S2_TOA collection by using consecutive filter methods on the collection. Before that, we will be pulling out a sentinel image from Sioux Falls, South Dakota. The tile number is T14TPP.

Lets, look for the day of 06/02/2019 because I know for sure that we have an image on that day

I also have a local copy of that image to tinker within QGIS. So, we can verify the results locally.

First, define some date ranges for our collection filter,

var start_date='2019-06-01';

var end_date='2019-06-03';


Then, apply the filters to the collection , then sort in ascending order and then just select the very first one from the sorted list

var s2_im=ee.Image(S2_TOA

                  .filterDate(start_date,end_date)

                  .filterBounds(point)

                  .sort('CLOUD_COVER')

                  .first());

 

Pack the visulaiztion parameters in a dictionary object

var trueColor={bands:['B4','B3','B2'],min:0,max:6000};

 

Use the addlayer method of Map class to add the image in the map

Map.addLayer(s2_im,trueColor,'True Color')

 

 print(s2_im,'all bands of the image')

//Map.addLayer(s2_im)

Clicking on the console tab we can explore bands available for the s2_im variable.

The expansion of the Image on the console yields the following. It has all the bands from Sentinel. 

Now, we will filter only the VNIR bands (bands which are Red (B4), Green (B3), Blue(B2) and NIR(B8) )

To select only VNIR bands we need to crate a list of the band IDs . After creating the list we need to pass the bandID  list to the ‘select’ method of the Image ‘s2_im’ object.

 

Now, to select only VNIR bands, make a variable that corresponds to the VNIR bands

var vnir_list=['B2','B3','B4','B8'];

var s2_im_vnir=s2_im.select(vnir_list);

var wavelengths=[0.490, 0.560, 0.665, 0.842];//optional for legend only

Check the s2_im  by printing to see if only vnir is selected

print(s2_im_vnir,'vnir filterd')

The output should be like the following after expansion in the console tab.

So, we have succesfully selected the only VNIR bands. The following part will do the tricks for extraction and plotting of the data--

Extarcting and plotting the VNIR data from the point geometry:

We need to configure the visualization option parameter first. To do that, we make a dictionary/object with the options.


var options={

title: 'VNIR raw data extarction from S2 using GEE (point Geometry)',

hAxis:{title:'Bands'},

vAxis:{title:'Reflectance(Raw DN)'},

lineWidth:2,

pointSize:7

};

 

 

This will create a chart in the console window when we print the chart object.

var chart=ui.Chart.image.regions(

  s2_im_vnir,point,null,null,null,vnir_list

  ).setOptions(options)

 print(chart,'Extract from point geometry')



Extarcting and plotting the VNIR data from the polygon geometry:

Similar method is used for reducing a polygon geometry. Also, we need to configure the visualization option parameter again here.

var options={

title: 'VNIR raw data extarction from S2 using GEE (Polygon Geometry)',

hAxis:{title:'Bands'},

vAxis:{title:'Reflectance(Raw DN)'},

lineWidth:2,

pointSize:7

};


var chart=ui.Chart.image.regions(

  s2_im_vnir,poly,ee.Reducer.mean(),null,null,vnir_list

  ).setOptions(options)

 

print(chart,'Extract from rectangle geometry')


The results are shown below:





Validation part: 
So the pixel values received from the google earth engine is plotted on the above picture.
To do the validation we will be using QGIS and a hand downloaded scene (same scene)
.Before that we want to export the point and polygon as kml/shp file. The process of doing that is fairly simple.

This snippet of comments are for extracting the geometry to kml or shape file. It will be exported to the google drive and can be  accessed through the account

First, make a collection of points.

var features = ee.FeatureCollection([

  ee.Feature(poly, {name: 'Point_of_interest_S2'})

]);

Export the FeatureCollection to a KML/SHP file.

Export.table.toDrive({

  collection: features,

  description:'rect_angle',

  fileFormat: 'SHP'

});

We can see these export code in the Tasks tab. It will not be exported to the google drive until we run the code from the Tasks list.

Opening and verification in QGIS:

I have already exported the point and polygon in the google drive and loaded in the local computer. So, I am doing a screenshot of  the QGIS window. and the QGIS window output is the following.


If we zoom in we can see the point and the polygon.


Now to validate the point geometry we need to zoom in to the point and get to the pixel level and get the DN value of the desired layer. To do that we need the 'value tool' for qgis. it can be installed as plugin if not available. We will need a 'zonal statistics'  plugin also for polygon verification. 

Doing a zoom to the pixel level and checking the DN value of that pixel through the value tool reveals that the Band2 or Blue Band value is 1643. Now, go back to the chart generated from GEE for point geometry,w here I intentionally put the cursor to on the B2(blue point) where it revealed 1643 also. So we extracted the right pixel.

Now, to verify the mean DN value of the polygon we need to run zonal statistics with the polygon shapefile and B2. To do that click on the zonal statistics toolbox from the 'processing toolbox section' and select the polygon layer and Blue band and process it.

Now, if we check the attribute table of the shape file we will see the mean value of all pixels inside the polygon for Blue band. 


The value we got in return is 1729.59810. Going back to the GEE returned calculation, we see the value is 1729.635 . The mean value is also matching 


The codes all together can also be found here.


No comments:

Post a Comment