Process satellite imagery

Need to add imagery to your map? This guide will teach you the basics of image processing for mapping, including an introduction to raster data and how to acquire, publish, and process raster imagery of our world.

Open source and at your fingertips. Let’s dive in.

The raster data type

Rasters are a pixel-based data format that efficiently represent continuous surfaces. The other common type of data you will come across in the world of mapping is vector data. Information in a raster is stored in a grid structure with each unit of information, or pixel, having the same size and shape but varying in value. Digital photographs are stored in this format, as are satellite images and orthophotography:

rasters-are-pixels

Raster formats are well suited to analyses that look at change over space and time. Each piece of information has an easily accessible location based on the grid. Using the location-based information that we have in the grid, we can access the same geographic location in two different georeferenced images and compare their values. By the end of this guide, you’ll have first-hand experience with this kind of image comparison.

Satellite rasters

When Earth-observing satellites take a picture, they read and record reflectance values collected from wavelengths along the electromagnetic spectrum.

The human eye can only see a small portion of the light-energy that constitutes the electromagnetic spectrum. This is called ‘visible light’, because our vision evolved to be most sensible where the Sun emits the most light, and broadly restricted to wavelengths that we make what we call red, green, and blue. Satellite sensors perceive a far wider range of the electromagnetic spectrum. The ability of sensors to collect information outside of our normal range of vision allows us to make visible the previously invisible.

The electromagnetic spectrum has a wide range, and it would be impractical for a sensor to collect information from all wavelengths at one time. Instead, different sensors prioritize the collection of information from different wavelengths of the spectrum. Each section of the spectrum that is captured and bucketed by a sensor is categorized as a band of information.

Bands of information, which vary in size, can be compiled into different types of composite images, each emphasizing a distinct physical property. In this guide, we will be creating true color composite images, prioritizing the natural look of land and water.

To make a true color composite, we’ll use the red, green and blue bands from images collected by Landsat 7; and we will do another composite with Landsat 8, to compare both. If you’re interested in learning more about different bands available for compositing, head over to our blog post on Putting Landsat 8 bands to work.

Getting started

To get started, you’ll need a few items. Full disclosure: this guide uses the command line.

Besides utilities like tar and brew, you’ll need the current versions of:

  • GDAL, a low-level GIS toolkit (brew install gdal if you don’t have it)

  • libgeotiff, to work with geotags (brew install libgeotiff if you don’t have it)

  • ImageMagick, an image processing package (brew install imagemagick -with--libtiff if you don’t have it)

  • Mapbox Studio Classic, our open-source map design studio

We’ll explain what’s going on in each step, so if you prefer other tools you’ll be able to translate the techniques.

Getting a scene

Satellite data is cut into scenes, which are near-square images covering an area that varies in size depending on the satellite and sensor that’s collecting the information. For Landsat images, which we’ll use in this guide, the coverage is about 170 × 185 km (105 × 115 mi) per scene. Think of a scene as a single frame from a camera.

There are several sources to find and download imagery:

If you’re a beginner, you should spend a few minutes with A Quick Guide to Earth Explorer for Landsat 8 from Robert Simmon, a data visualization expert at NASA’s wonderful Earth Observatory.

Setting the scene

For this guide, we want to compare historical and present-day images of Dubai in the United Arab Emirates. Dubai experienced an incredible amount of growth in the early 2000s. Let’s see if we can find an image before and after its economic boom.

Downloading the scenes

Following Simmon’s guide, we’ve found that Landsat 7 scene LE71600432000140SGS00 from late May 2000 and Landsat 8 scene LC81600432014154LGN00 from early June 2014 are cloud-free images that suit our needs.

Landsat Archive Landsat Scene Identifier
L7 ETM+ SLC-on (1999-2003) LE71600432000140SGS00
L8 OLI/TIRS LC81600432014154LGN00

Let’s download the Level 1 product of each of these scenes from EarthExplorer:

  1. From EarthExplorer, click the Data Sets tab on the left-side panel.
  2. Click Landsat Archive to reveal more options.
  3. Check the Landsat’s respective archive name (shown in the table above).
  4. Click Additional Criteria.
  5. Enter the Landsat Scene Identifier (shown in the table above).
  6. Click Results.
  7. Click the download icon (green arrow pointing down).
  8. Click Download Level 1 Product. (Note that Landsat 7 requires processing; click the What’s This? link for further instructions. It may take up to 24 hours to process.)

In the bundle

The Level 1 products come as .tar.gz files ranging in size from 700 megabytes to a gigabyte for Landsat 8 and 200-400 megabytes for Landsat 7. The compressed file will unpack into a directory of 13 items, mostly TIFF images, each with an unwieldy-looking name starting with the 21-digit scene ID. This is the bundle.

Here’s what you need to know about the bundle:

  • The images with names ending in digits are the data for those bands. For example, LC80120542013154LGN00_B9.TIF is Band 9’s readout.
  • The data is Level 1T terrain corrected, meaning that it’s been filtered to account for some sensor variations and for distortions caused by hills and valleys. The correction is not perfect, for example because the elevation data that it’s referenced against is a little coarse, but L1T is a good starting point for most uses.
  • Each image is aligned with the others, so you know that a given pixel at x, y in one band’s data represents exactly the same point in space as the corresponding pixel in another image in the same bundle. The exception is Landsat 8’s Band 8, the pan band, which is at twice the resolution. (The alignment does not carry between different bundles that share path/row numbers. They can vary by 10 km or so.)

These properties tell us enough to dive right into creating images.

Processing for true color

Now that we have the images we want, let’s get down to processing. We’ll work entirely from the command line for this section. For now, pull up a Terminal window and navigate to the folder where your Landsat 7 bundle lives.

Landsat 7

On Landsat 7 the red, green, and blue bands are the ones numbered 3, 2, and 1. Let’s first create a variable for our scene ID, and then make reprojected versions into Web Mercator (EPSG:3857) of each of the bands we will be working with:

#establish an id variable
id="LE71600432000140SGS00"

#use a for loop to reproject each of the bands you will be working with.
for BAND in {3,2,1}; do
gdalwarp -t_srs EPSG:3857 $id"_B"$BAND.tif $BAND-projected.tif;done

Now let’s merge them into an RGB image with gdal_merge.py:

#merge the three reprojected band images
#into a single composite image
gdal_merge.py -v -separate -of GTiff -co PHOTOMETRIC=RGB -o $id-RGB.tif 3-projected.tif 2-projected.tif 1-projected.tif

Now that we have our RGB image, LE71600432000140SGS00-RGB.tif, let’s take a look at it:

rasters-L7-dubai

It’s a little hazy, but overall the color, saturation, and contrast are pretty good. Let’s use ImageMagick to tweak them just a bit:

#modulate default brightness, increased saturation.
#sigmoidal contrast increased contrast, default middle brightness
convert -channel RGB -modulate 100,150 -sigmoidal-contrast 4x50% $id-RGB.tif $id-RGB-cc.tif

Incidentally, don’t worry when convert prints about a dozen warnings like this one:

convert: Unknown field with tag 33550 (0x830e) encountered.
`TIFFReadDirectory' @ warning/tiff.c/TIFFWarnings/824.

ImageMagick has a number of functions, only a few of which we are using here:

  • First we use the -channel flag to make sure we are working with all three bands at the same time. You may also work with single bands at a time, which we will do with our next image.
  • The -modulate flag allows us to adjust overall brightness and saturation.
  • The -sigmoidal-contrast flag is for contrast and establishing how far into your value-range you’d like middle brightness to be.

Let’s see how our image looks now:

rasters-L7-dubai-cc

Those small corrections make for a much better looking image, but we’re not finished yet.

First, we will build overviews with gdaladdo. Building in overviews is an important optimization step when working with images that are many megabytes or gigabytes in size. The gdaladdo process inserts reduced resolution image versions inside the original image file, which are used in place of the full resolution data at low zoom levels. Building overviews helps maintain a quick load time by not requiring more detail to load than is necessary at each zoom level.

We also need to account for the fact that ImageMagick isn’t geo-aware. The geodata gets knocked off of your file when you process it with ImageMagick commands. To re-attach the geodata let’s pull a copy of it from one of our images that has not be processed with ImageMagick. In the example, we’ll use the version of band 3 that we reprojected with gdalwarp just a few steps back. It doesn’t matter which image you choose to pull your geodata from, so long as its projection is consistent with the output you desire.

#use a cubic downsampling method to add overviews
#(other interpolation methods are available)
gdaladdo -r cubic $id-RGB-cc.tif 2 4 8 10 12

#create a TIFF worldfile from the requested image
listgeo -tfw 3-projected.tif

The listgeo command takes the geodata from your input image and rewrites it to a file with the same basename but the suffix .tfw (“TIFF worldfile”), which we’ll rename to match the file that we’re re-georeferencing:

#rename the requested TIFF worldfile to have same name
#as the file that needs to be re-georeferenced
mv 3-projected.tfw $id-RGB-cc.tfw

GDAL knows to look for a matching .tfw if it sees a TIFF that isn’t internally georeferenced:

#apply spatial reference system to non-georeferenced file
#(GDAL knows to look for files with matching names)
gdal_edit.py -a_srs EPSG:3857 $id-RGB-cc.tif

For the last step, let’s remove the black background.

gdalwarp -srcnodata 0 -dstalpha $id-RGB-cc.tif $id-RGB-cc-2.tif

rasters-L7-dubai-cc

Landsat 7 is all set. Let’s work on Landsat 8 next.

Landsat 8

We will repeat much of the process we just went through for our Landsat 8 image with one important addition. We need to convert the image to 8-bit. Landsat 8 data is collected by a sensor that has a higher radiometric resolution than the sensor on Landsat 7. This means that Landsat 8 data comes in a 16-bit format, while Landsat 7 comes as 8-bit. Our styling tool, Studio, will only render images up to the 8-bit resolution, so we need to convert our 16-bit Landsat 8 bands into 8-bit.

To reformat the image, we need to take the value range in the 16-bit image (0 - 65,535) and squish it down to the same value range that is found in an 8-bit image (0 - 255). The next steps follow the same process as with our Landsat 7 image, but we’ll apply gdal_translate commands to rescale the bands to the 8-bit format.

Here’s the full script:

#establish an id variable
id="LC81600432014154LGN00"

#use a for loop to reproject each of the bands you will be working with.
for BAND in {4,3,2}; do
 gdalwarp -t_srs EPSG:3857 $id"_B"$BAND.tif $BAND-projected.tif;done

#translate each of your bands into the 8-bit format with default settings of -ot and -scale
gdal_translate -ot Byte -scale 0 65535 0 255 4-projected{,-scaled}.tif
gdal_translate -ot Byte -scale 0 65535 0 255 3-projected{,-scaled}.tif
gdal_translate -ot Byte -scale 0 65535 0 255 2-projected{,-scaled}.tif

#merge the three reprojected band images into a single composite image
gdal_merge.py -v -ot Byte -separate -of GTiff -co PHOTOMETRIC=RGB -o $id-RGB-scaled.tif 4-projected-scaled.tif 3-projected-scaled.tif 2-projected-scaled.tif

#color corrections in blue bands to deal with haze factor,
#and across all brands for brightness, contrast and saturation
convert -channel B -gamma 1.05 -channel RGB -sigmoidal-contrast 20,20% -modulate 100,150 $id-RGB-scaled.tif $id-RGB-scaled-cc.tif

#use a cubic downsampling method to add overview
#(other interpolation methods are available)
gdaladdo -r cubic $id-RGB-scaled-cc.tif 2 4 8 10 12

#call the TIFF worldfile for the requested image,
#change name of file to match file needing georeference,
#and apply georeference
listgeo -tfw 3-projected.tif
mv 3-projected.tfw $id-RGB-scaled-cc.tfw
gdal_edit.py -a_srs EPSG:3857 $id-RGB-scaled-cc.tif

#remove black background
gdalwarp -srcnodata 0 -dstalpha $id-RGB-scaled-cc.tif $id-RGB-scaled-cc-2.tif

This long script prepares a single Landsat 8 image doing the warp, merge, and color adjustment.

And here are the original RGB and the color corrected output:

rasters-L8-dubai rasters-L8-dubai-cc

Working in Mapbox Studio Classic

Importing

Now, let’s import our images into Mapbox Studio Classic. Repeat the following steps for both images since each image will have its own data source.

  1. Open up Mapbox Studio Classic and start a new project.
  2. Select Blank source under New source.
  3. Click New layer.
  4. Browse and select Landsat image.
  5. Click Done.
  6. In Settings update the name of the datasource to reflect the landsat number.
  7. Click Upload to Mapbox.

Alternatively, you can upload the images directly on Mapbox.com. Check out Carol’s post: Raster imagery to Mapbox in one step.

You can find the uploaded sources on your Data page. Each source will have its own map ID. You can use the map ID to reference the data when using developer tools, such as Mapbox.js. This is great for layering your data on top of a map.

Styling

If you’d like to make adjustments to your image or layer it on a map, you can style it.

From the data source’s Settings panel, click Create style from source and then click Remote source. You’ll see the image on a blank map with a few CartoCSS styles to the right.

Let’s edit the CartoCSS to give the image the best rendering:

#_image {
  raster-opacity: 1;
  raster-scaling: lanczos;
}

Let’s also remove the Map {...} styles to keep the background transparent.

In the Docs panel, you can find more CartoCSS properties and values to style the image.

Once you’ve finished styling, publish your map by clicking Upload to Mapbox in the Settings panel. Your published map will appear in your Styles page with its map ID.

Making it interactive

For the final step, we’ll add some interaction to our map tiles with the swipe layers function.

You can copy and paste the code snippet from this example and replace the map IDs with your own two unique map IDs for the Dubai images we’ve just created.

Next, center the map on Dubai. This can be done with the map.setView function. We used LatLng 25.0859, 55.1720 and zoom level 11. Feel free to play around with it to find an initial view that you think is effective in engaging the user.

Your final product will be a compelling view of landscape change that you can embed in any webpage you’d like:

Mission complete

Nice work! You just acquired your own raster data, processed it for best visual impact, created a compelling display structure, and published your product to the web.

If you enjoyed this guide then you should check out Georeferencing imagery next. We also recommend that you browse our Mapbox.js examples for ideas to push your imagery even further.

Additional questions? Ask our support team or learn more about How Mapbox Works.