Extracting the central strip from LANDSAT 7 imagery

Extracting the central strip from LANDSAT 7 imagery

February 8, 2016

Here is a simple Python code to extract the central strip from Landsat 7 imagery (SLC-off),  that is not affected by the SLC failure. The algorithm shrinks the striping zones through a morphological filter (erosion) and creates a new shapefile AOI that extracts the desired raster extent without striping effects. The code is based on Python for ArcGIS (arcpy) – so you require a ArcGIS license.

General steps:

  1. Loop through all Landsat 7 data folders
  2. Stack bands for each image
  3. Create a mask
  4. Erode the mask by 20 pixels
  5. Convert the mask to polygon
  6. Create a minimum bounding box
  7. Clip the original raster through the bbox

 

import arcpy
from arcpy.sa import *

import sys,os

#  Environment settings (Activate Spatial Analyst, Overwrite Outputs allowed and TIFF compression is LZW)
arcpy.CheckOutExtension("spatial")
arcpy.env.overwriteOutput = True
arcpy.env.compression = 'LZW'

# this is your main directory with all unzipped Landsat datasets
 rootdir = "D:\\DATA\\Landsat7\\"

# create scratch folder "temp" 
temp = "D:\\DATA\\temp\\"

# loop through directory with all unzipped Landsat 7 folders
 for files in os.listdir(rootdir):   
    files = os.path.join(rootdir, files)   
    
    # for each loop the subdir "files" is now the current workspace 
    # (e.g. LE71520322015157-SC20160224113319) that contains the Landsat bands
    arcpy.env.workspace = files  
    rasters = arcpy.ListRasters("*", "TIF")  
    
    # create empty array
    stack_liste = []  
    # loop through all rasters in subdir
    for raster in rasters:   

        image = arcpy.Raster(raster) 
        name  = image.name 
        index = name.split("_")[0]  

        # fill up the array only with the actual spectral bands        
        sr = "_sr_band"  
        if sr in raster:   
            stack_liste.append(raster)             

    # now stack all bands within the array
    stack_name = files + "\\" + index + "_stack.tif"    
    arcpy.CompositeBands_management(stack_liste, stack_name)  

    # convert the image stack to a mask by logical operation with an absurd value that will result in an output "0"
    con = EqualTo(stack_name, 123456789)  

    # now shrink the raster mask with value "0" by 20 pixels
    shrink = temp + "shrink"  
    shrinking = Shrink(con, 20, 0) 
    shrinking.save(shrink)  

    zone = temp + "zone.shp" 
    bbox = temp + "bbox.shp"  

    # conver the shrunk mask to polygon and create a minimum bounding box
    arcpy.RasterToPolygon_conversion(shrink, zone, "NO_SIMPLIFY", "VALUE") 
    arcpy.MinimumBoundingGeometry_management(zone, bbox, "RECTANGLE_BY_WIDTH", "NONE")  

    # now use that bounding box as a mask to cut out the central nadir strip from the original stack
    # Final result 
    extract = files + "\\" + index + "_aoi.tif"  
    ExtractByMask = arcpy.sa.ExtractByMask(stack_name, bbox) 
    ExtractByMask.save(extract)

 

follow us and share it on:

you may also like:

Starkregen in Bayern: Beobachtungen und Dokumentation zählen

Starkregen in Bayern: Beobachtungen und Dokumentation zählen

Starkregenereignisse treten immer häufiger lokal, kurzfristig und mit hoher Intensität auf. Innerhalb weniger Stunden können sie erhebliche Überschwemmungen und Schäden verursachen. Um solche Ereignisse künftig besser zu verstehen und die wissenschaftliche Grundlage...

Seeing the World in Points: Lidar Course for the EAGLEs

Seeing the World in Points: Lidar Course for the EAGLEs

Lidar has a funny way of sneaking up on you. You think you know what it is, a laser that measures distance, fine, but then someone shows you a point cloud of a forest canopy with individual branches floating in 3D space and suddenly you realize there's a whole...

RTL covers EORC: TV Crew Films MONID Habitrack Fieldwork

RTL covers EORC: TV Crew Films MONID Habitrack Fieldwork

A bit of extra excitement at EORC recently, an RTL television crew showed up to film a segment on the MONID Habitrack project, and Dr. Ariane Droin was right in the middle of it, walking them through what Earth Observation actually brings to the table for a project...

Ticks from Above: UAS Fieldwork for the MONID Habitrack Project

Ticks from Above: UAS Fieldwork for the MONID Habitrack Project

Forest edges are tricky places. They're where woodland meets open ground, where light and shade trade off every few meters, and where, it turns out, ticks tend to do really well. That last bit is exactly why Dr. Ariane Droin, Sofica Garcia de Leon, Dr. Jakob...

Course on urban EO by Michael Wurm

Course on urban EO by Michael Wurm

Walk through any city and you pick up on things that are hard to put a number on. The noise of a main road, the heat that sits between buildings in summer, the question of whether that little park around the corner is really enough green space for the whole...

EireR R package: unified gateway to Irish geospatial data

EireR R package: unified gateway to Irish geospatial data

Anyone who's tried to do geospatial work across the whole island of Ireland knows the headache. Ireland is one island geographically, but it's split across two jurisdictions, the Republic and Northern Ireland, and each one runs its own data infrastructure. Different...

Share This