Getting started

Add the crate

[dependencies]
rwarp = "0.3"

For unreleased changes on main, depend on the repository:

[dependencies]
rwarp = { git = "https://github.com/hypertidy/rwarp" }

Choose a CRS backend

cargo test                                            # libproj
cargo test --no-default-features --features proj4rs   # pure Rust
cargo build --no-default-features --features proj4rs --target wasm32-unknown-unknown

Or bring your own: GenImgProjTransformer::with_backend takes any pair of CrsTransform implementations (one each way), for example a JavaScript PROJ build reached from WebAssembly.

A first warp

This plans and warps one 256 x 256 tile in a local Lambert azimuthal equal-area projection from a Web Mercator tile pyramid. It is examples/tile_laea.rs; run it with cargo run --example tile_laea.

use rwarp::source_window::compute_source_window;
use rwarp::transform::GenImgProjTransformer;
use rwarp::warp::{warp, Nodata, ResampleAlg, SourceBuffer, WarpOptions};

const WEBMERC_HALF: f64 = 20037508.342789244;
const LAEA_TAS: &str = "+proj=laea +lat_0=-42 +lon_0=147 +datum=WGS84 +units=m +no_defs";

fn main() -> Result<(), String> {
    // Source: Web Mercator at zoom 8, as one virtual raster.
    let n = 256 * (1 << 8);
    let px = 2.0 * WEBMERC_HALF / n as f64;
    let src_gt = [-WEBMERC_HALF, px, 0.0, WEBMERC_HALF, 0.0, -px];
    let src_size = [n, n];

    // Target: a 256 x 256 tile, 400 km across, centred on Tasmania.
    let extent = 400_000.0;
    let dst_gt = [-extent / 2.0, extent / 256.0, 0.0, extent / 2.0, 0.0, -extent / 256.0];

    // 1. Transformer: target pixel -> source pixel, approximated as gdalwarp does.
    let exact = GenImgProjTransformer::new("EPSG:3857", src_gt, LAEA_TAS, dst_gt)?;
    let t = ApproxTransformer::new(exact, 0.125);

    // 2. Plan: which source pixels does the tile need? No pixels are read.
    let win = compute_source_window(&t, [0, 0], [256, 256], src_size, 1)
        .ok_or("tile does not map onto the source")?;
    println!(
        "source window: offset ({}, {}), size {} x {} pixels",
        win.xoff, win.yoff, win.xsize, win.ysize
    );

    // 3. Read: fetch that window (here, synthesise it).
    let (w, h) = (win.xsize as usize, win.ysize as usize);
    let data: Vec<f32> = (0..w * h).map(|i| (win.xoff as usize + i % w) as f32).collect();
    let src = SourceBuffer {
        data: &data,
        ncol: w,
        nrow: h,
        nbands: 1,
        col_off: win.xoff as usize,
        row_off: win.yoff as usize,
    };

    // 4. Warp.
    let opts = WarpOptions::new(ResampleAlg::Bilinear)
        .with_src_nodata(Nodata::Value(f64::NAN));
    let out = warp(&t, &src, 256, 256, &opts);
    let valid = out.iter().filter(|v| !v.is_nan()).count();
    println!("warped {} of {} target pixels", valid, out.len());
    Ok(())
}

The four steps are the whole pipeline:

1. Transformer. Maps target pixels to source pixels. Wrapping it in ApproxTransformer with a 0.125 pixel error, as gdalwarp does by default, transforms a few points per scanline exactly and interpolates the rest.
2. Plan. compute_source_window finds the source pixels the tile needs, without reading any.
3. Read. Fetch that window however your data is stored.
4. Warp. warp resamples into the target, with options for the method, kernel scaling, nodata and more (see Resampling).