martin_tile_utils/
lib.rs

1#![doc = include_str!("../README.md")]
2
3// This code was partially adapted from https://github.com/maplibre/mbtileserver-rs
4// project originally written by Kaveh Karimi and licensed under MIT/Apache-2.0
5
6use std::f64::consts::PI;
7use std::fmt::{Display, Formatter, Result};
8
9/// circumference of the earth in meters
10pub const EARTH_CIRCUMFERENCE: f64 = 40_075_016.685_578_5;
11/// circumference of the earth in degrees
12pub const EARTH_CIRCUMFERENCE_DEGREES: u32 = 360;
13
14/// radius of the earth in meters
15pub const EARTH_RADIUS: f64 = EARTH_CIRCUMFERENCE / 2.0 / PI;
16
17pub const MAX_ZOOM: u8 = 30;
18
19mod decoders;
20pub use decoders::*;
21mod rectangle;
22pub use rectangle::{TileRect, append_rect};
23
24#[derive(Debug, Copy, Clone, Hash, PartialEq, Eq)]
25pub struct TileCoord {
26    pub z: u8,
27    pub x: u32,
28    pub y: u32,
29}
30
31pub type TileData = Vec<u8>;
32pub type Tile = (TileCoord, Option<TileData>);
33
34impl Display for TileCoord {
35    fn fmt(&self, f: &mut Formatter<'_>) -> Result {
36        if f.alternate() {
37            write!(f, "{}/{}/{}", self.z, self.x, self.y)
38        } else {
39            write!(f, "{},{},{}", self.z, self.x, self.y)
40        }
41    }
42}
43
44impl TileCoord {
45    /// Checks provided coordinates for validity
46    /// before constructing [`TileCoord`] instance.
47    ///
48    /// Check [`Self::new_unchecked`] if you are sure that your inputs are possible.
49    #[must_use]
50    pub fn new_checked(z: u8, x: u32, y: u32) -> Option<TileCoord> {
51        Self::is_possible_on_zoom_level(z, x, y).then_some(Self { z, x, y })
52    }
53
54    /// Constructs [`TileCoord`] instance from arguments without checking that the tiles can exist.
55    ///
56    /// Check [`Self::new_checked`] if you are unsure if your inputs are possible.
57    #[must_use]
58    pub fn new_unchecked(z: u8, x: u32, y: u32) -> TileCoord {
59        Self { z, x, y }
60    }
61
62    /// Checks that zoom `z` is plausibily small and `x`/`y` is possible on said zoom level
63    #[must_use]
64    pub fn is_possible_on_zoom_level(z: u8, x: u32, y: u32) -> bool {
65        if z > MAX_ZOOM {
66            return false;
67        }
68
69        let side_len = 1_u32 << z;
70        x < side_len && y < side_len
71    }
72}
73
74#[derive(Clone, Copy, Debug, PartialEq, Eq)]
75pub enum Format {
76    Gif,
77    Jpeg,
78    Json,
79    Mvt,
80    Png,
81    Webp,
82}
83
84impl Format {
85    #[must_use]
86    pub fn parse(value: &str) -> Option<Self> {
87        Some(match value.to_ascii_lowercase().as_str() {
88            "gif" => Self::Gif,
89            "jpg" | "jpeg" => Self::Jpeg,
90            "json" => Self::Json,
91            "pbf" | "mvt" => Self::Mvt,
92            "png" => Self::Png,
93            "webp" => Self::Webp,
94            _ => None?,
95        })
96    }
97
98    /// Get the `format` value as it should be stored in the `MBTiles` metadata table
99    #[must_use]
100    pub fn metadata_format_value(self) -> &'static str {
101        match self {
102            Self::Gif => "gif",
103            Self::Jpeg => "jpeg",
104            Self::Json => "json",
105            // QGIS uses `pbf` instead of `mvt` for some reason
106            Self::Mvt => "pbf",
107            Self::Png => "png",
108            Self::Webp => "webp",
109        }
110    }
111
112    #[must_use]
113    pub fn content_type(&self) -> &str {
114        match *self {
115            Self::Gif => "image/gif",
116            Self::Jpeg => "image/jpeg",
117            Self::Json => "application/json",
118            Self::Mvt => "application/x-protobuf",
119            Self::Png => "image/png",
120            Self::Webp => "image/webp",
121        }
122    }
123
124    #[must_use]
125    pub fn is_detectable(self) -> bool {
126        match self {
127            Self::Png | Self::Jpeg | Self::Gif | Self::Webp => true,
128            // TODO: Json can be detected, but currently we only detect it
129            //       when it's not compressed, so to avoid a warning, keeping it as false for now.
130            //       Once we can detect it inside a compressed data, change it to true.
131            Self::Mvt | Self::Json => false,
132        }
133    }
134}
135
136impl Display for Format {
137    fn fmt(&self, f: &mut Formatter<'_>) -> Result {
138        f.write_str(match *self {
139            Self::Gif => "gif",
140            Self::Jpeg => "jpeg",
141            Self::Json => "json",
142            Self::Mvt => "mvt",
143            Self::Png => "png",
144            Self::Webp => "webp",
145        })
146    }
147}
148
149#[derive(Clone, Copy, Debug, Hash, PartialEq, Eq)]
150pub enum Encoding {
151    /// Data is not compressed, but it can be
152    Uncompressed = 0b0000_0000,
153    /// Some formats like JPEG and PNG are already compressed
154    Internal = 0b0000_0001,
155    Gzip = 0b0000_0010,
156    Zlib = 0b0000_0100,
157    Brotli = 0b0000_1000,
158    Zstd = 0b0001_0000,
159}
160
161impl Encoding {
162    #[must_use]
163    pub fn parse(value: &str) -> Option<Self> {
164        Some(match value.to_ascii_lowercase().as_str() {
165            "none" => Self::Uncompressed,
166            "gzip" => Self::Gzip,
167            "zlib" => Self::Zlib,
168            "brotli" => Self::Brotli,
169            "zstd" => Self::Zstd,
170            _ => None?,
171        })
172    }
173
174    #[must_use]
175    pub fn content_encoding(&self) -> Option<&str> {
176        match *self {
177            Self::Uncompressed | Self::Internal => None,
178            Self::Gzip => Some("gzip"),
179            Self::Zlib => Some("deflate"),
180            Self::Brotli => Some("br"),
181            Self::Zstd => Some("zstd"),
182        }
183    }
184
185    #[must_use]
186    pub fn is_encoded(self) -> bool {
187        match self {
188            Self::Uncompressed | Self::Internal => false,
189            Self::Gzip | Self::Zlib | Self::Brotli | Self::Zstd => true,
190        }
191    }
192}
193
194#[derive(Clone, Copy, Debug, PartialEq, Eq)]
195pub struct TileInfo {
196    pub format: Format,
197    pub encoding: Encoding,
198}
199
200impl TileInfo {
201    #[must_use]
202    pub fn new(format: Format, encoding: Encoding) -> Self {
203        Self { format, encoding }
204    }
205
206    /// Try to figure out the format and encoding of the raw tile data
207    #[must_use]
208    #[allow(clippy::enum_glob_use)]
209    pub fn detect(value: &[u8]) -> Option<Self> {
210        use Encoding::*;
211        use Format::*;
212
213        // TODO: Make detection slower but more accurate:
214        //  - uncompress gzip/zlib/... and run detection again. If detection fails, assume MVT
215        //  - detect json inside a compressed data
216        //  - json should be fully parsed
217        //  - possibly keep the current `detect()` available as a fast path for those who may need it
218        Some(match value {
219            // Compressed prefixes assume MVT content
220            v if v.starts_with(b"\x1f\x8b") => Self::new(Mvt, Gzip),
221            v if v.starts_with(b"\x78\x9c") => Self::new(Mvt, Zlib),
222            v if v.starts_with(b"\x89\x50\x4E\x47\x0D\x0A\x1A\x0A") => Self::new(Png, Internal),
223            v if v.starts_with(b"\x47\x49\x46\x38\x39\x61") => Self::new(Gif, Internal),
224            v if v.starts_with(b"\xFF\xD8\xFF") => Self::new(Jpeg, Internal),
225            v if v.starts_with(b"RIFF") && v.len() > 8 && v[8..].starts_with(b"WEBP") => {
226                Self::new(Webp, Internal)
227            }
228            v if v.starts_with(b"{") => Self::new(Json, Uncompressed),
229            _ => None?,
230        })
231    }
232
233    #[must_use]
234    pub fn encoding(self, encoding: Encoding) -> Self {
235        Self { encoding, ..self }
236    }
237}
238
239impl From<Format> for TileInfo {
240    fn from(format: Format) -> Self {
241        Self::new(
242            format,
243            match format {
244                Format::Png | Format::Jpeg | Format::Webp | Format::Gif => Encoding::Internal,
245                Format::Mvt | Format::Json => Encoding::Uncompressed,
246            },
247        )
248    }
249}
250
251impl Display for TileInfo {
252    fn fmt(&self, f: &mut Formatter<'_>) -> Result {
253        write!(f, "{}", self.format.content_type())?;
254        if let Some(encoding) = self.encoding.content_encoding() {
255            write!(f, "; encoding={encoding}")?;
256        } else if self.encoding != Encoding::Uncompressed {
257            f.write_str("; uncompressed")?;
258        }
259        Ok(())
260    }
261}
262
263/// Convert longitude and latitude to a tile (x,y) coordinates for a given zoom
264#[must_use]
265#[allow(clippy::cast_possible_truncation)]
266#[allow(clippy::cast_sign_loss)]
267pub fn tile_index(lng: f64, lat: f64, zoom: u8) -> (u32, u32) {
268    let tile_size = EARTH_CIRCUMFERENCE / f64::from(1_u32 << zoom);
269    let (x, y) = wgs84_to_webmercator(lng, lat);
270    let col = (((x - (EARTH_CIRCUMFERENCE * -0.5)).abs() / tile_size) as u32).min((1 << zoom) - 1);
271    let row = ((((EARTH_CIRCUMFERENCE * 0.5) - y).abs() / tile_size) as u32).min((1 << zoom) - 1);
272    (col, row)
273}
274
275/// Convert min/max XYZ tile coordinates to a bounding box values.
276///
277/// The result is `[min_lng, min_lat, max_lng, max_lat]`
278///
279/// # Panics
280/// Panics if `zoom` is greater than [`MAX_ZOOM`].
281#[must_use]
282pub fn xyz_to_bbox(zoom: u8, min_x: u32, min_y: u32, max_x: u32, max_y: u32) -> [f64; 4] {
283    assert!(zoom <= MAX_ZOOM, "zoom {zoom} must be <= {MAX_ZOOM}");
284
285    let tile_length = EARTH_CIRCUMFERENCE / f64::from(1_u32 << zoom);
286
287    let left_down_bbox = tile_bbox(min_x, max_y, tile_length);
288    let right_top_bbox = tile_bbox(max_x, min_y, tile_length);
289
290    let (min_lng, min_lat) = webmercator_to_wgs84(left_down_bbox[0], left_down_bbox[1]);
291    let (max_lng, max_lat) = webmercator_to_wgs84(right_top_bbox[2], right_top_bbox[3]);
292    [min_lng, min_lat, max_lng, max_lat]
293}
294
295#[allow(clippy::cast_lossless)]
296fn tile_bbox(x: u32, y: u32, tile_length: f64) -> [f64; 4] {
297    let min_x = EARTH_CIRCUMFERENCE * -0.5 + x as f64 * tile_length;
298    let max_y = EARTH_CIRCUMFERENCE * 0.5 - y as f64 * tile_length;
299
300    [min_x, max_y - tile_length, min_x + tile_length, max_y]
301}
302
303/// Convert bounding box to a tile box `(min_x, min_y, max_x, max_y)` for a given zoom
304#[must_use]
305pub fn bbox_to_xyz(left: f64, bottom: f64, right: f64, top: f64, zoom: u8) -> (u32, u32, u32, u32) {
306    let (min_col, min_row) = tile_index(left, top, zoom);
307    let (max_col, max_row) = tile_index(right, bottom, zoom);
308    (min_col, min_row, max_col, max_row)
309}
310
311/// Compute precision of a zoom level, i.e. how many decimal digits of the longitude and latitude are relevant
312#[must_use]
313#[allow(clippy::cast_possible_truncation, clippy::cast_sign_loss)]
314pub fn get_zoom_precision(zoom: u8) -> usize {
315    assert!(zoom < MAX_ZOOM, "zoom {zoom} must be <= {MAX_ZOOM}");
316    let lng_delta = webmercator_to_wgs84(EARTH_CIRCUMFERENCE / f64::from(1_u32 << zoom), 0.0).0;
317    let log = lng_delta.log10() - 0.5;
318    if log > 0.0 { 0 } else { -log.ceil() as usize }
319}
320
321/// transform [`WebMercator`](https://epsg.io/3857) to [WGS84](https://epsg.io/4326)
322// from https://github.com/Esri/arcgis-osm-editor/blob/e4b9905c264aa22f8eeb657efd52b12cdebea69a/src/OSMWeb10_1/Utils/WebMercator.cs
323#[must_use]
324pub fn webmercator_to_wgs84(x: f64, y: f64) -> (f64, f64) {
325    let lng = (x / EARTH_RADIUS).to_degrees();
326    let lat = f64::atan(f64::sinh(y / EARTH_RADIUS)).to_degrees();
327    (lng, lat)
328}
329
330/// transform [WGS84](https://epsg.io/4326) to [`WebMercator`](https://epsg.io/3857)
331// from https://github.com/Esri/arcgis-osm-editor/blob/e4b9905c264aa22f8eeb657efd52b12cdebea69a/src/OSMWeb10_1/Utils/WebMercator.cs
332#[must_use]
333pub fn wgs84_to_webmercator(lon: f64, lat: f64) -> (f64, f64) {
334    let x = lon * PI / 180.0 * EARTH_RADIUS;
335
336    let y_sin = lat.to_radians().sin();
337    let y = EARTH_RADIUS / 2.0 * ((1.0 + y_sin) / (1.0 - y_sin)).ln();
338
339    (x, y)
340}
341
342#[cfg(test)]
343mod tests {
344    #![allow(clippy::unreadable_literal)]
345
346    use std::fs::read;
347
348    use Encoding::{Internal, Uncompressed};
349    use Format::{Jpeg, Json, Png, Webp};
350    use approx::assert_relative_eq;
351    use rstest::rstest;
352
353    use super::*;
354
355    fn detect(path: &str) -> Option<TileInfo> {
356        TileInfo::detect(&read(path).unwrap())
357    }
358
359    #[allow(clippy::unnecessary_wraps)]
360    fn info(format: Format, encoding: Encoding) -> Option<TileInfo> {
361        Some(TileInfo::new(format, encoding))
362    }
363
364    #[test]
365    fn test_data_format_png() {
366        assert_eq!(detect("./fixtures/world.png"), info(Png, Internal));
367    }
368
369    #[test]
370    fn test_data_format_jpg() {
371        assert_eq!(detect("./fixtures/world.jpg"), info(Jpeg, Internal));
372    }
373
374    #[test]
375    fn test_data_format_webp() {
376        assert_eq!(detect("./fixtures/dc.webp"), info(Webp, Internal));
377        assert_eq!(TileInfo::detect(br"RIFF"), None);
378    }
379
380    #[test]
381    fn test_data_format_json() {
382        assert_eq!(
383            TileInfo::detect(br#"{"foo":"bar"}"#),
384            info(Json, Uncompressed)
385        );
386    }
387
388    #[rstest]
389    #[case(-180.0, 85.0511, 0, (0,0))]
390    #[case(-180.0, 85.0511, 1, (0,0))]
391    #[case(-180.0, 85.0511, 2, (0,0))]
392    #[case(0.0, 0.0, 0, (0,0))]
393    #[case(0.0, 0.0, 1, (1,1))]
394    #[case(0.0, 0.0, 2, (2,2))]
395    #[case(0.0, 1.0, 0, (0,0))]
396    #[case(0.0, 1.0, 1, (1,0))]
397    #[case(0.0, 1.0, 2, (2,1))]
398    fn test_tile_colrow(
399        #[case] lng: f64,
400        #[case] lat: f64,
401        #[case] zoom: u8,
402        #[case] expected: (u32, u32),
403    ) {
404        assert_eq!(
405            expected,
406            tile_index(lng, lat, zoom),
407            "{lng},{lat}@z{zoom} should be {expected:?}"
408        );
409    }
410
411    #[rstest]
412    // you could easily get test cases from maptiler: https://www.maptiler.com/google-maps-coordinates-tile-bounds-projection/#4/-118.82/71.02
413    #[case(0, 0, 0, 0, 0, [-180.0,-85.0511287798066,180.0,85.0511287798066])]
414    #[case(1, 0, 0, 0, 0, [-180.0,0.0,0.0,85.0511287798066])]
415    #[case(5, 1, 1, 2, 2, [-168.75,81.09321385260837,-146.25,83.97925949886205])]
416    #[case(5, 1, 3, 2, 5, [-168.75,74.01954331150226,-146.25,81.09321385260837])]
417    fn test_xyz_to_bbox(
418        #[case] zoom: u8,
419        #[case] min_x: u32,
420        #[case] min_y: u32,
421        #[case] max_x: u32,
422        #[case] max_y: u32,
423        #[case] expected: [f64; 4],
424    ) {
425        let bbox = xyz_to_bbox(zoom, min_x, min_y, max_x, max_y);
426        assert_relative_eq!(bbox[0], expected[0], epsilon = f64::EPSILON * 2.0);
427        assert_relative_eq!(bbox[1], expected[1], epsilon = f64::EPSILON * 2.0);
428        assert_relative_eq!(bbox[2], expected[2], epsilon = f64::EPSILON * 2.0);
429        assert_relative_eq!(bbox[3], expected[3], epsilon = f64::EPSILON * 2.0);
430    }
431
432    #[rstest]
433    #[case(0, (0, 0, 0, 0))]
434    #[case(1, (0, 1, 0, 1))]
435    #[case(2, (0, 3, 0, 3))]
436    #[case(3, (0, 7, 0, 7))]
437    #[case(4, (0, 14, 1, 15))]
438    #[case(5, (0, 29, 2, 31))]
439    #[case(6, (0, 58, 5, 63))]
440    #[case(7, (0, 116, 11, 126))]
441    #[case(8, (0, 233, 23, 253))]
442    #[case(9, (0, 466, 47, 507))]
443    #[case(10, (1, 933, 94, 1014))]
444    #[case(11, (3, 1866, 188, 2029))]
445    #[case(12, (6, 3732, 377, 4059))]
446    #[case(13, (12, 7465, 755, 8119))]
447    #[case(14, (25, 14931, 1510, 16239))]
448    #[case(15, (51, 29863, 3020, 32479))]
449    #[case(16, (102, 59727, 6041, 64958))]
450    #[case(17, (204, 119455, 12083, 129917))]
451    #[case(18, (409, 238911, 24166, 259834))]
452    #[case(19, (819, 477823, 48332, 519669))]
453    #[case(20, (1638, 955647, 96665, 1039339))]
454    #[case(21, (3276, 1911295, 193331, 2078678))]
455    #[case(22, (6553, 3822590, 386662, 4157356))]
456    #[case(23, (13107, 7645181, 773324, 8314713))]
457    #[case(24, (26214, 15290363, 1546649, 16629427))]
458    #[case(25, (52428, 30580726, 3093299, 33258855))]
459    #[case(26, (104857, 61161453, 6186598, 66517711))]
460    #[case(27, (209715, 122322907, 12373196, 133035423))]
461    #[case(28, (419430, 244645814, 24746393, 266070846))]
462    #[case(29, (838860, 489291628, 49492787, 532141692))]
463    #[case(30, (1677721, 978583256, 98985574, 1064283385))]
464    fn test_box_to_xyz(#[case] zoom: u8, #[case] expected_xyz: (u32, u32, u32, u32)) {
465        let actual_xyz = bbox_to_xyz(
466            -179.43749999999955,
467            -84.76987877980656,
468            -146.8124999999996,
469            -81.37446385260833,
470            zoom,
471        );
472        assert_eq!(
473            actual_xyz, expected_xyz,
474            "zoom {zoom} does not have te right xyz"
475        );
476    }
477
478    #[rstest]
479    // test data via https://epsg.io/transform#s_srs=4326&t_srs=3857
480    #[case((0.0,0.0), (0.0,0.0))]
481    #[case((30.0,0.0), (3339584.723798207,0.0))]
482    #[case((-30.0,0.0), (-3339584.723798207,0.0))]
483    #[case((0.0,30.0), (0.0,3503549.8435043753))]
484    #[case((0.0,-30.0), (0.0,-3503549.8435043753))]
485    #[case((38.897957,-77.036560), (4330100.766138651, -13872207.775755845))] // white house
486    #[case((-180.0,-85.0), (-20037508.342789244, -19971868.880408566))]
487    #[case((180.0,85.0), (20037508.342789244, 19971868.880408566))]
488    #[case((0.026949458523585632,0.08084834874097367), (3000.0, 9000.0))]
489    fn test_coordinate_syste_conversion(
490        #[case] wgs84: (f64, f64),
491        #[case] webmercator: (f64, f64),
492    ) {
493        // epsg produces the expected values with f32 precision, grrr..
494        let epsilon = f64::from(f32::EPSILON);
495
496        let actual_wgs84 = webmercator_to_wgs84(webmercator.0, webmercator.1);
497        assert_relative_eq!(actual_wgs84.0, wgs84.0, epsilon = epsilon);
498        assert_relative_eq!(actual_wgs84.1, wgs84.1, epsilon = epsilon);
499
500        let actual_webmercator = wgs84_to_webmercator(wgs84.0, wgs84.1);
501        assert_relative_eq!(actual_webmercator.0, webmercator.0, epsilon = epsilon);
502        assert_relative_eq!(actual_webmercator.1, webmercator.1, epsilon = epsilon);
503    }
504
505    #[rstest]
506    #[case(0..11, 0)]
507    #[case(11..14, 1)]
508    #[case(14..17, 2)]
509    #[case(17..21, 3)]
510    #[case(21..24, 4)]
511    #[case(24..27, 5)]
512    #[case(27..30, 6)]
513    fn test_get_zoom_precision(
514        #[case] zoom: std::ops::Range<u8>,
515        #[case] expected_precision: usize,
516    ) {
517        for z in zoom {
518            let actual_precision = get_zoom_precision(z);
519            assert_eq!(
520                actual_precision, expected_precision,
521                "Zoom level {z} should have precision {expected_precision}, but was {actual_precision}"
522            );
523        }
524    }
525
526    #[test]
527    fn test_tile_coord_zoom_range() {
528        for z in 0..=MAX_ZOOM {
529            assert!(TileCoord::is_possible_on_zoom_level(z, 0, 0));
530            assert_eq!(
531                TileCoord::new_checked(z, 0, 0),
532                Some(TileCoord { z, x: 0, y: 0 })
533            );
534        }
535        assert!(!TileCoord::is_possible_on_zoom_level(MAX_ZOOM + 1, 0, 0));
536        assert_eq!(TileCoord::new_checked(MAX_ZOOM + 1, 0, 0), None);
537    }
538
539    #[test]
540    fn test_tile_coord_new_checked_xy_for_zoom() {
541        assert!(TileCoord::is_possible_on_zoom_level(5, 0, 0));
542        assert_eq!(
543            TileCoord::new_checked(5, 0, 0),
544            Some(TileCoord { z: 5, x: 0, y: 0 })
545        );
546        assert!(TileCoord::is_possible_on_zoom_level(5, 31, 31));
547        assert_eq!(
548            TileCoord::new_checked(5, 31, 31),
549            Some(TileCoord { z: 5, x: 31, y: 31 })
550        );
551        assert!(!TileCoord::is_possible_on_zoom_level(5, 31, 32));
552        assert_eq!(TileCoord::new_checked(5, 31, 32), None);
553        assert!(!TileCoord::is_possible_on_zoom_level(5, 32, 31));
554        assert_eq!(TileCoord::new_checked(5, 32, 31), None);
555    }
556
557    #[test]
558    /// Any (u8, u32, u32) values can be put inside [`TileCoord`], of course, but some
559    /// functions may panic at runtime (e.g. [`mbtiles::invert_y_value`]) if they are impossible,
560    /// so let's not do that.
561    fn test_tile_coord_new_unchecked() {
562        assert_eq!(
563            TileCoord::new_unchecked(u8::MAX, u32::MAX, u32::MAX),
564            TileCoord {
565                z: u8::MAX,
566                x: u32::MAX,
567                y: u32::MAX
568            }
569        );
570    }
571
572    #[test]
573    fn xyz_format() {
574        let xyz = TileCoord { z: 1, x: 2, y: 3 };
575        assert_eq!(format!("{xyz}"), "1,2,3");
576        assert_eq!(format!("{xyz:#}"), "1/2/3");
577    }
578}
martin_tile_utils/lib.rs

martin_tile_utils/
lib.rs
