1.Extract_HLS_L30 .R

Syntax is easy. Knowing which libraries to use – not so much!

Here’s a quick description of each Python library used for this program:

• pandas: Provides data structures such as DataFrames. Aliased as pd .

• rasterio: Provides functions for processing geospatial raster data.

• geopandas: Extends pandas to support spatial data. Aliased as gpd.

• os: Provides operating system features such as changing directories.

• glob: Provides functions for listing files that match specific patterns.

• numpy: Provides functions for multi-dimensional arrays and matrices. Aliased as np.

• rasterstats: Summarizes geospatial raster datasets based on vector geometries. We only need the zonal_stats function from this library.

Defines a function that takes the sensor (L30 or S30) and NIR_band ( B05 or B8A) as parameters.

Changes to the base directory and gets a list of the tiles.

Creates a dataframe with the number of rows in fields but with no columns.

Note the simple syntax in Python. The single colon : means “all rows.” 0:0 means “no columns.”

Creates a loop to iterate through the tiles.

Note that the Python code could be written as for tile in tiles. However, the i variable gives us the opportunity of showing the status of the program, as described later.

Changes to the directory of the current tile.
Note that the Python enumerate function allows us to grab both the tile and its index into the tiles list.
Therefore, we don’t need to refer to tiles[i], although that would work just as well.

Grabs the files that match the pattern HLS.{sensor}.

Python’s f strings work similarly to R’s glue function.

Creates a loop to iterate through each tile's scenes.

Shows the current processing status.

We need to account for the fact that Python counts from 0 instead of 1, so it outputs j + 1 and i + 1 instead of j and i.

Changes to the scene/date folder.

Grabs the filenames of the NIR and B04 bands.

Note how Python’s glob.glob function simplifies this.

Loads the raster stack of the B04 and NIR bands.

Note that rast.stk is a valid variable name in R. However, the dot (.) makes it a syntax error in Python.

Per PEP 8, variable names “should be lowercase, with words separated by underscores.”

Sets the names of the raster stack based on file names.

The R code creates references to cloud_mask. So cloud_mask0, cloud_mask64, and cloud_mask128 essentially become “aliases” of cloud_mask.

This is an optimization technique, which prevents R from copying data that might not change.

However, the data does change in the following code. When that happens, R will automatically change the references to copies of cloud_mask.

The Python code makes this explicit by creating copies of cloud_mask to begin with.

Sets values not matching the conditions to NaN.

Merges the cloud mask layers for 0, 64, and 128.

Applies the cloud mask to the raster stack

CONTINUE FROM HERE

Reads the buffered shapefile.

Plots the shapefile.

edgecolor isn’t required in the Python code. It simply makes the output consistent with the R output.