import numpy as np import rasterio from rasterio.transform import from_origin from rasterio.warp import transform_bounds from pipeline.transform import noise_overlay_tiles from pipeline.transform.noise_overlay_tiles import RasterInfo, _read_noise_tile def _write_corridor_raster(path, nodata=-96.0): """A small EPSG:27700 raster: a column of 70 dB cells adjacent to genuine 0.0 (quiet) cells. Bilinear blending of the 0 cells would fabricate a halo of intermediate dB values between 0 and 70.""" # 8x8 grid: leftmost two columns are 70 dB, the rest are genuine quiet 0.0. data = np.zeros((8, 8), dtype=np.float32) data[:, 0:2] = 70.0 # Place one true nodata cell to make sure it is also masked out. data[0, 7] = nodata # 10m cells anchored somewhere inside England's BNG extent. left = 300_000.0 top = 300_080.0 transform = from_origin(left, top, 10.0, 10.0) with rasterio.open( path, "w", driver="GTiff", height=data.shape[0], width=data.shape[1], count=1, dtype=data.dtype, crs="EPSG:27700", transform=transform, nodata=nodata, ) as dataset: dataset.write(data, 1) return path def test_read_noise_tile_does_not_fabricate_halo(tmp_path): raster_path = _write_corridor_raster(tmp_path / "corridor.tif") with rasterio.open(raster_path) as dataset: bounds_27700 = dataset.bounds bounds_mercator = transform_bounds( dataset.crs, noise_overlay_tiles.WEB_MERCATOR_CRS, *bounds_27700, densify_pts=21, ) info = RasterInfo(path=raster_path, bounds_mercator=bounds_mercator) # Render at high resolution so any bilinear halo would surface as # intermediate dB values along the corridor/quiet seam. tile_size = 64 tile = _read_noise_tile([info], bounds_mercator, tile_size) finite = tile[np.isfinite(tile)] # Every finite cell must be the genuine corridor value (~70). There must be # NO fabricated halo strictly between 0 and 70. halo = finite[(finite > 0.0) & (finite < 70.0 - 1e-3)] assert halo.size == 0, f"fabricated halo values present: {np.unique(halo)}" # Sanity: the corridor itself must still be rendered. assert finite.size > 0 assert np.all(finite >= 70.0 - 1e-3) def test_read_noise_tile_preserves_peak_under_downsample(tmp_path): # 8x8 EPSG:27700 raster: a single loud 75 dB cell in a 50 dB field. # Downsampling into a smaller tile with bilinear would dilute the peak # (arithmetic dB averaging); Resampling.max must keep the worst-case dB. data = np.full((8, 8), 50.0, dtype=np.float32) data[4, 4] = 75.0 transform = from_origin(300_000.0, 300_080.0, 10.0, 10.0) raster_path = tmp_path / "peak.tif" with rasterio.open( raster_path, "w", driver="GTiff", height=data.shape[0], width=data.shape[1], count=1, dtype=data.dtype, crs="EPSG:27700", transform=transform, nodata=-96.0, ) as dataset: dataset.write(data, 1) with rasterio.open(raster_path) as dataset: bounds_mercator = transform_bounds( dataset.crs, noise_overlay_tiles.WEB_MERCATOR_CRS, *dataset.bounds, densify_pts=21, ) info = RasterInfo(path=raster_path, bounds_mercator=bounds_mercator) # Render the 8x8 source into a 4x4 tile: this downsamples, so bilinear # would average the 75 dB peak away. tile = _read_noise_tile([info], bounds_mercator, 4) finite = tile[np.isfinite(tile)] assert finite.size > 0 # The loud peak must survive the downsample (max, not arithmetic mean). assert finite.max() >= 75.0 - 1e-3, f"peak diluted to {finite.max()}" # Max resampling must never invent a value louder than the source. assert finite.max() <= 75.0 + 1e-3