ArcGIS Python API Fire Analysis

[KrishV]

Hi Everyone,

I'm trying to use California wildfires 2017 - Thomas Fire analysis Sample provided in ArcGIS Python API Samples. I have used sentinel imagery to do the same analysis but I always end up getting black imagery (All 0s). Can you please suggest where I'm doing wrong?

import arcgis
from arcgis import *
from arcgis.mapping import MapImageLayer

gis = GIS(profile="krishna_dev") 

#Sentinel2 
sentinel_item = gis.content.search('title:Sentinel2', 'Imagery Layer', outside_org=True)[0]
sentinel = sentinel_item.layers[0]
sentinel_item

aoi = {'spatialReference': {'latestWkid': 3857, 'wkid': 102100}, 'type': 'extent', 
       'xmax': 9028482.969, 'xmin': 8908018.212, 'ymax': 2478700.211, 'ymin': 2248166.133}

arcgis.env.analysis_extent = {'xmax': 9028482.969, 'xmin': 8908018.212, 'ymax': 2478700.211, 'ymin': 2248166.133,
                              "spatialReference":{"wkid":102100,"latestWkid":3857}}
sentinel.extent = aoi

import pandas as pd
from datetime import datetime


selected = sentinel.filter_by(where="(Category = 1)",
                             time=[datetime(2020, 4, 3), datetime(2020, 4, 15)],
geometry=arcgis.geometry.filters.intersects(aoi))

df = selected.query(out_fields="acquisitiondate, name, cloudcover", 
                    order_by_fields="acquisitiondate").sdf
df['acquisitiondate'] = pd.to_datetime(df['acquisitiondate'], unit='ms')
df.tail()

prefire = sentinel.filter_by('OBJECTID in ( 10127407, 10187367)') # 2020-04-03
prefire = prefire.blend()

midfire = sentinel.filter_by('OBJECTID in (10304749, 10304746)') # 2020-04-13 
midfire = midfire.blend()

from arcgis.raster.functions import *

extracted_band = extract_band(midfire, [13,12,4])
extracted_band_prefire = extract_band(prefire, [13,12,4])

nbr_prefire  = band_arithmetic(prefire, "(b9 - b13) / (b9 + b13) + 2000")
nbr_postfire = band_arithmetic(midfire, "(b9 - b13) / (b9 + b13) + 2000")

nbr_diff = nbr_prefire - nbr_postfir

nbr_diff

burnt_areas = colormap(remap(nbr_diff,
                             input_ranges=[0.1,  0.27,  # low severity 
                                           0.27, 0.44,  # medium severity
                                           0.44, 0.66,  # moderate severity
                                           0.66, 1.00], # high severity burn
                             output_values=[1, 2, 3, 4],                    
                             no_data_ranges=[-1, 0.1], astype='u8'), 
                             colormap=[[4, 0xFF, 0xC3, 0], [3, 0xFA, 0x8E, 0], [2, 0xF2, 0x55, 0], [1, 0xE6, 0,    0]])

burnt_areas

Thanks,

Krish

[MehdiPira1]

Hi @KrishV 

There is no Band 13 in Sentinel2 and your NBR formula

"(b9 - b13) / (b9 + b13) + 2000"

is different from that of esri.

The correct formula for NBR is 

"(b5 - b7) / (b5 + b7+1000)"

 Besides that, you have used 

time=[datetime(2020, 4, 3), datetime(2020, 4, 15)]

this period of time which means there might be no fire occurred between these dates.

If you visualize prefire and midfire, you might be able to check if a fire has occurred in the images or not prior to further calculations.

Cheers

Mehdi

[KrishV]

Hi @MehdiPira1 ,

Did you get a chance to check on this?

Thanks,

Krish

[KrishV]

Hi @MehdiPira1 ,

I have taken the formula from below site:

https://developers.arcgis.com/python/sample-notebooks/wildfire-analysis-using-sentinel-2-imagery/

I have used below formula also but I don't get any desire results:

"(b5 - b7) / (b5 + b7+1000)"

I have also visualized the midfire image and I see the burnt areas

Thanks,

Krish

[MehdiPira1]

Hi @KrishV ,

That's bizarre!

Did you try to apply that manually in ArcGIS Pro to see if the process gives the same result or not?

[KrishV · Accepted solution]

Hi @MehdiPira1 ,

I was getting result in ArcGIS Pro using Sentinel Images but same thing is not happening using ArcGIS Python API if I use Sentinel Imagery available in Living Atlas.

However, I have managed to extract burnt areas using Landsat Imagery available in Living Atlas

Thanks,

Krish