How fast is my city sprawling — in built-up land and at night?

A growing city shows up from space two ways. The ground turns from green or bare dirt into grey rooftops and roads, and the night gets brighter as new streets, homes, and shops switch their lights on. Nighttime lights are exactly that: a satellite measurement of how much light shines up from the Earth’s surface after dark, in a unit called radiance (basically “brightness per pixel”). Hold up one clear night against another a few years later and you can measure how much brighter your city has become.

Which data to use, and how

Use this dataset: Black Marble nighttime lights (short name VNP46A2), NASA’s daily map of how bright the Earth’s surface is at night, at about ~500 m per pixel, going back to 2012. The handy part: it’s already gap-filled, so pixels hidden by clouds on a given night are patched in and you don’t start with a Swiss-cheese image. The field you want is Gap_Filled_DNB_BRDF-Corrected_NTL, in units of nW/cm²/sr (a standard brightness unit; bigger number means brighter).

Then do four steps. The runnable code further down does all of this; this is just the idea.

1. Pick your city: a small latitude/longitude box around it — and two dates a few years apart.
2. Crop and average: pull the brightness for every pixel inside your box on one clear night and take the average, giving one brightness number for that night.
3. Repeat for the later year: do the same on a night in the same season a few years later (same season matters; see the gotchas).
4. Compare: subtract, or take the percent change. Up means your city got brighter (usually growth). You can also map where it brightened most to spot the expanding edge.

That’s the whole method. Below, “How a scientist answers this” names the exact tools, and the Run it cell does it live on synthetic data so you can see each step work.

What you can find out

• How bright your city is on any clear night, back to 2012. The gap-filled field already patches cloud-blocked pixels for you.
• How much brighter it got over a few years. Difference two clear nights and map where the brightness rose the most (often the city’s expanding edge).
• Where new lit land appeared. Count the share of pixels above a “lit” brightness threshold each year and watch that fraction climb as the urban footprint spreads.
• How nighttime growth lines up with built-up land. Overlay daytime built-up cover from HLSL30 or the free GHSL settlement layer to separate “just brighter” from “genuinely newly built”.
• Verified locally: over the Bengaluru box (77.4–77.8 °E, 12.8–13.2 °N), averaging five clear January nights, mean nighttime brightness rose from about 28.3 nW/cm²/sr in 2018 to 33.4 in 2024, roughly an 18% increase (a single best-night comparison gave +23%). That box was already ~96–100% “lit” in 2018, so the growth showed up as existing areas getting brighter rather than dark land switching on. Read lit-area and brightness together.

What it can’t tell you

• What the light is for. The satellite sees brightness, not whether it’s homes, factories, highways, or one new stadium. You can’t read land use from radiance alone; you’d need other maps.
• Population or economic growth directly. Lights correlate with both, but a brighter pixel isn’t a headcount or a GDP figure. Treat it as a rough proxy, never a measurement.
• Whether new land is actually built up. A brighter night could just be more (or more efficient) lighting on the same footprint. To confirm real construction you need daytime imagery (HLSL30, GHSL).
• Fine within-block detail. At ~500 m, a single pixel blends a whole city block together. You can’t isolate one building or one street.
• A reading through clouds or a bright moon. Even after gap-filling, persistently cloudy stretches and full-moon nights add noise. That’s why you always compare clear nights under similar moonlight.

Gotchas to watch for

• Compare the same season. Festival lighting (Diwali, Christmas) and winter haze can make a night look much brighter without any real growth, and that’s the easiest way to manufacture a “trend” that isn’t there. Pick the same month in both years: pit January against January, never January against a festival week.
• Skip cloudy and bright-moon nights. Gap-filling patches cloud holes, but a sky that’s cloudy for weeks, or lit by a full moon, still adds noise to the brightness number. Average several clear nights per year and try to match the moon phase.
• Read brightness and lit-area together. Track both the average brightness and the share of pixels above a “lit” threshold. A city that’s already fully lit (like the Bengaluru example) grows by getting brighter, not by adding dark land, so a “lit pixel count” alone would miss it.
• Use the right field and crop first. VNP46A2 ships as daily HDF5 tiles (a NASA file format), and in the current (v002) files the data sits under a grid called VIIRS_Grid_DNB_2d, with each tile carrying its own latitude/longitude arrays. Skip the crop and you’ll average a huge area and get a meaningless number, so use the tile’s own lat/lon to cut down to your box before you compute.

The real-data code

The Run it cell above runs this method on synthetic data with no login. Below is the whole recipe against the real archive: search Black Marble tiles, open the HDF5, crop to your box, average several clear January nights per year, compare brightness and lit-area, and plot the change. It needs a free Earthdata Login (one free account, no fee).

import os, re, warnings, earthaccess, h5py, numpy as np
import matplotlib.pyplot as plt
warnings.filterwarnings("ignore")

# load Earthdata creds from .env without `source` (passwords can break the shell)
for line in open(".env"):
    m = re.match(r'\s*(?:export\s+)?([A-Z0-9_]+)\s*=\s*(.*)\s*$', line)
    if m: os.environ.setdefault(m.group(1), m.group(2).strip().strip('"').strip("'"))
earthaccess.login(strategy="environment")   # free Earthdata Login

W, S, E, N = 77.4, 12.8, 77.8, 13.2          # your city (Bengaluru) as W,S,E,N
GRID = "HDFEOS/GRIDS/VIIRS_Grid_DNB_2d/Data Fields"
VAR  = GRID + "/Gap_Filled_DNB_BRDF-Corrected_NTL"   # the gap-filled brightness field
LIT  = 5.0                                    # nW/cm2/sr threshold counting a pixel as "lit"

def crop_one_night(start, end):
    # 1-2. find Black Marble tiles for this date + box, download one, open the HDF5, crop to box
    g = earthaccess.search_data(short_name="VNP46A2",
                                temporal=(start, end), bounding_box=(W, S, E, N))
    if not g:
        return None
    fn = str(earthaccess.download(g[:1], local_path="/tmp/ntl")[0])
    with h5py.File(fn, "r") as f:
        lat, lon = f[GRID + "/lat"][:], f[GRID + "/lon"][:]
        yi = np.where((lat >= S) & (lat <= N))[0]   # rows inside your box
        xi = np.where((lon >= W) & (lon <= E))[0]   # columns inside your box
        ds  = f[VAR]
        arr = ds[yi.min():yi.max()+1, xi.min():xi.max()+1].astype("float32")
        arr[arr == float(ds.attrs["_FillValue"][0])] = np.nan   # drop "no data" pixels
    return arr

def year_stats(year):
    # 3. average several clear January nights so one cloudy/bright-moon night can't skew it
    nights = [crop_one_night(f"{year}-01-{d:02d}", f"{year}-01-{d:02d}")
              for d in (8, 10, 12, 14, 16)]
    nights = [a for a in nights if a is not None]
    stack  = np.dstack(nights)                       # rows x cols x nights
    mean_scene = np.nanmean(stack, axis=2)           # per-pixel mean across clear nights
    brightness = float(np.nanmean(mean_scene))       # one brightness number for the year
    lit_frac   = float(np.nanmean(mean_scene > LIT)) # share of box that is "lit"
    return mean_scene, brightness, lit_frac

# 4. compare two years a few seasons-matched Januaries apart
scene18, b18, lit18 = year_stats(2018)
scene24, b24, lit24 = year_stats(2024)
pct = (b24 - b18) / b18 * 100

verdict = "brighter (likely growth)" if pct > 0 else "dimmer"
print(f"Bengaluru nighttime brightness {b18:.1f} -> {b24:.1f} nW/cm2/sr "
      f"({pct:+.0f}%), {verdict}; lit area {lit18*100:.0f}% -> {lit24*100:.0f}%")

# map WHERE it brightened most — the expanding edge
fig, (ax1, ax2) = plt.subplots(1, 2, figsize=(9, 4))
ax1.imshow(scene24, cmap="inferno", vmin=0, vmax=60); ax1.set_title("2024 brightness")
d = ax2.imshow(scene24 - scene18, cmap="RdBu_r", vmin=-30, vmax=30)
ax2.set_title("change 2018->2024"); fig.colorbar(d, ax=ax2, label="d nW/cm2/sr")
for ax in (ax1, ax2): ax.set_xticks([]); ax.set_yticks([])
plt.tight_layout(); plt.show()

# Before you trust it: a near-fully-lit box grows by getting *brighter*, not by adding dark land,
# so read brightness and lit-area together, and confirm new construction with daytime HLSL30 / GHSL.

Where the data comes from

The nighttime lights come from NASA’s VIIRS instrument (flying on the Suomi-NPP satellite), processed into the Black Marble product and distributed by the LAADS archive, reachable through one free Earthdata Login. The optional daytime “is it really built up?” check uses HLSL30 (also NASA/Earthdata) or GHSL, a free built-up settlement layer from the European Union. Separate, but free.

Sources

• VNP46A2 Black Marble user guide: https://ladsweb.modaps.eosdis.nasa.gov/missions-and-measurements/products/VNP46A2/
• Black Marble project: https://blackmarble.gsfc.nasa.gov/
• GHSL (EU built-up settlement layer): https://ghsl.jrc.ec.europa.eu/