Skip to content

Commit cb1106e

Browse files
kylebarronkylebarron
authored
feat(geotiff): generate TileMatrixSet from GeoTIFF instance (#235)
### Change list - Add new `generateTileMatrixSet` API to create a `TileMatrixSet` object from a GeoTIFF's internal layout. - Add test to mimic Python test in https://github.com/developmentseed/async-geotiff/blob/e835e9531049a2631d326478c002ef3d8c0296f3/tests/test_tms.py#L15-L53 - Remove `GeoTIFF.epsg`. All CRS access should go through `GeoTIFF.crs` - Update readme of morecantile package
1 parent 20b3ace commit cb1106e

File tree

12 files changed

+354
-29
lines changed

12 files changed

+354
-29
lines changed

‎packages/affine/src/affine.ts‎

Lines changed: 30 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -83,3 +83,33 @@ export function invert([sa, sb, sc, sd, se, sf]: Affine): Affine {
8383

8484
return [ra, rb, -sc * ra - sf * rb, rd, re, -sc * rd - sf * re];
8585
}
86+
87+
/** Get the 'a' component of an Affine transform. */
88+
export function a(affine: Affine): number {
89+
return affine[0];
90+
}
91+
92+
/** Get the 'b' component of an Affine transform. */
93+
export function b(affine: Affine): number {
94+
return affine[1];
95+
}
96+
97+
/** Get the 'c' component of an Affine transform. */
98+
export function c(affine: Affine): number {
99+
return affine[2];
100+
}
101+
102+
/** Get the 'd' component of an Affine transform. */
103+
export function d(affine: Affine): number {
104+
return affine[3];
105+
}
106+
107+
/** Get the 'e' component of an Affine transform. */
108+
export function e(affine: Affine): number {
109+
return affine[4];
110+
}
111+
112+
/** Get the 'f' component of an Affine transform. */
113+
export function f(affine: Affine): number {
114+
return affine[5];
115+
}

‎packages/geotiff/package.json‎

Lines changed: 3 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -47,6 +47,8 @@
4747
},
4848
"dependencies": {
4949
"@cogeotiff/core": "^9.1.2",
50-
"@developmentseed/affine": "workspace:^"
50+
"@developmentseed/affine": "workspace:^",
51+
"@developmentseed/morecantile": "workspace:^",
52+
"uuid": "^13.0.0"
5153
}
5254
}

‎packages/geotiff/src/geotiff.ts‎

Lines changed: 0 additions & 10 deletions
Original file line numberDiff line numberDiff line change
@@ -147,16 +147,6 @@ export class GeoTIFF {
147147
return this._crs;
148148
}
149149

150-
/** EPSG code from GeoTIFF tags, or null if not set.
151-
*
152-
* See also {@link GeoTIFF.crs} for the full PROJJSON definition, which should
153-
* always be available, even when an EPSG code is not explicitly stored in the
154-
* tags.
155-
*/
156-
get epsg(): number | null {
157-
return this.image.epsg;
158-
}
159-
160150
/** Image width in pixels. */
161151
get width(): number {
162152
return this.image.size.width;

‎packages/geotiff/src/index.ts‎

Lines changed: 1 addition & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -5,3 +5,4 @@ export { decode, registry } from "./decode/api.js";
55
export { GeoTIFF } from "./geotiff.js";
66
export { Overview } from "./overview.js";
77
export type { Tile } from "./tile.js";
8+
export { generateTileMatrixSet } from "./tile-matrix-set.js";

‎packages/geotiff/src/tile-matrix-set.ts‎

Lines changed: 166 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,166 @@
1+
import type { Affine } from "@developmentseed/affine";
2+
import * as affine from "@developmentseed/affine";
3+
import type {
4+
BoundingBox,
5+
CRS,
6+
TileMatrix,
7+
TileMatrixSet,
8+
} from "@developmentseed/morecantile";
9+
import { metersPerUnit } from "@developmentseed/morecantile";
10+
import { v4 as uuidv4 } from "uuid";
11+
import type { ProjJson } from "./crs.js";
12+
import type { GeoTIFF } from "./geotiff.js";
13+
14+
/**
15+
* A minimal projection definition compatible with what wkt-parser returns.
16+
*
17+
* This type extracts only the partial properties we need from the full
18+
* wkt-parser output.
19+
*/
20+
interface ProjectionDefinition {
21+
datum?: {
22+
/** Semi-major axis of the ellipsoid. */
23+
a: number;
24+
};
25+
a?: number;
26+
to_meter?: number;
27+
units?: string;
28+
}
29+
30+
const SCREEN_PIXEL_SIZE = 0.28e-3;
31+
32+
function buildCrs(crs: number | ProjJson): CRS {
33+
if (typeof crs === "number") {
34+
return {
35+
uri: `http://www.opengis.net/def/crs/EPSG/0/${crs}`,
36+
};
37+
}
38+
39+
// @ts-expect-error - typing issues between different projjson definitions.
40+
return {
41+
wkt: crs,
42+
};
43+
}
44+
45+
/**
46+
* Build a TileMatrix entry for a single resolution level.
47+
*/
48+
function buildTileMatrix(
49+
id: string,
50+
transform: Affine,
51+
mpu: number,
52+
cornerOfOrigin: "bottomLeft" | "topLeft",
53+
tileWidth: number,
54+
tileHeight: number,
55+
width: number,
56+
height: number,
57+
): TileMatrix {
58+
return {
59+
id,
60+
scaleDenominator: (affine.a(transform) * mpu) / SCREEN_PIXEL_SIZE,
61+
cellSize: affine.a(transform),
62+
cornerOfOrigin,
63+
pointOfOrigin: [affine.c(transform), affine.f(transform)],
64+
tileWidth,
65+
tileHeight,
66+
matrixWidth: Math.ceil(width / tileWidth),
67+
matrixHeight: Math.ceil(height / tileHeight),
68+
};
69+
}
70+
71+
/**
72+
* Generate a Tile Matrix Set from a GeoTIFF file.
73+
*
74+
* Produces one TileMatrix per overview (coarsest first) plus a final entry
75+
* for the full-resolution level. The GeoTIFF must be tiled.
76+
*
77+
* This requires a crs definition that includes a `units` property, so that we
78+
* can convert pixel sizes to physical screen units. Use [`wkt-parser`] to parse
79+
* a WKT string or PROJJSON object, then pass the result as the `crs` argument.
80+
*
81+
* [`wkt-parser`]: https://github.com/proj4js/wkt-parser
82+
*
83+
* @see https://docs.ogc.org/is/17-083r4/17-083r4.html
84+
*/
85+
export function generateTileMatrixSet(
86+
geotiff: GeoTIFF,
87+
crs: ProjectionDefinition,
88+
{ id = uuidv4() }: { id?: string } = {},
89+
): TileMatrixSet {
90+
const bbox = geotiff.bbox;
91+
const tr = geotiff.transform;
92+
93+
// Perhaps we should allow metersPerUnit to take any string
94+
const crsUnit = crs.units as
95+
| "m"
96+
| "metre"
97+
| "meter"
98+
| "meters"
99+
| "foot"
100+
| "us survey foot"
101+
| "degree"
102+
| undefined;
103+
104+
if (!crsUnit) {
105+
throw new Error(`CRS definition must include "units" property`);
106+
}
107+
108+
const semiMajorAxis = crs.a || crs.datum?.a;
109+
const mpu = metersPerUnit(crsUnit, { semiMajorAxis });
110+
const cornerOfOrigin: "bottomLeft" | "topLeft" =
111+
affine.e(tr) > 0 ? "bottomLeft" : "topLeft";
112+
113+
const tileMatrices: TileMatrix[] = [];
114+
115+
// Overviews are sorted finest-to-coarsest; reverse to emit coarsest first.
116+
const overviewsCoarseFirst = [...geotiff.overviews].reverse();
117+
118+
for (let idx = 0; idx < overviewsCoarseFirst.length; idx++) {
119+
const overview = overviewsCoarseFirst[idx]!;
120+
tileMatrices.push(
121+
buildTileMatrix(
122+
String(idx),
123+
overview.transform,
124+
mpu,
125+
cornerOfOrigin,
126+
overview.tileWidth,
127+
overview.tileHeight,
128+
overview.width,
129+
overview.height,
130+
),
131+
);
132+
}
133+
134+
// Full-resolution level is appended last.
135+
if (!geotiff.isTiled) {
136+
throw new Error("GeoTIFF must be tiled to generate a TMS.");
137+
}
138+
139+
tileMatrices.push(
140+
buildTileMatrix(
141+
String(geotiff.overviews.length),
142+
tr,
143+
mpu,
144+
cornerOfOrigin,
145+
geotiff.tileWidth,
146+
geotiff.tileHeight,
147+
geotiff.width,
148+
geotiff.height,
149+
),
150+
);
151+
152+
const tmsCrs = buildCrs(geotiff.crs);
153+
const boundingBox: BoundingBox = {
154+
lowerLeft: [bbox[0], bbox[1]],
155+
upperRight: [bbox[2], bbox[3]],
156+
crs: tmsCrs,
157+
};
158+
159+
return {
160+
title: "Generated TMS",
161+
id,
162+
crs: tmsCrs,
163+
boundingBox,
164+
tileMatrices,
165+
};
166+
}

‎packages/geotiff/tests/crs.test.ts‎

Lines changed: 1 addition & 4 deletions
Original file line numberDiff line numberDiff line change
@@ -3,16 +3,13 @@ import wktParser from "wkt-parser";
33
import { loadGeoTIFF } from "./helpers.js";
44

55
describe("test CRS", () => {
6-
it("can fetch EPSG CRS from epsg.io", async () => {
6+
it("returns EPSG code", async () => {
77
const geotiff = await loadGeoTIFF(
88
"uint8_rgb_deflate_block64_cog",
99
"rasterio",
1010
);
1111
const crs = geotiff.crs;
1212
expect(crs).toEqual(4326);
13-
14-
const epsg = geotiff.epsg;
15-
expect(epsg).toBe(4326);
1613
});
1714
});
1815

‎packages/geotiff/tests/geotiff.test.ts‎

Lines changed: 0 additions & 13 deletions
Original file line numberDiff line numberDiff line change
@@ -178,19 +178,6 @@ describe("GeoTIFF", () => {
178178
expect(geo.nodata).toBe(-9999);
179179
});
180180

181-
it("exposes epsg", async () => {
182-
const primary = mockImage({
183-
width: 100,
184-
height: 100,
185-
origin: [0, 0, 0],
186-
resolution: [1, -1, 0],
187-
epsg: 4326,
188-
});
189-
const tiff = mockTiff([primary]);
190-
const geo = await GeoTIFF.fromTiff(tiff);
191-
expect(geo.epsg).toBe(4326);
192-
});
193-
194181
it("exposes bbox", async () => {
195182
const primary = mockImage({
196183
width: 100,

‎packages/geotiff/tests/tile-matrix-set.test.ts‎

Lines changed: 128 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,128 @@
1+
import { describe, expect, it } from "vitest";
2+
import wktParser from "wkt-parser";
3+
import { generateTileMatrixSet } from "../src/tile-matrix-set.js";
4+
import { loadGeoTIFF } from "./helpers.js";
5+
6+
const EPSG_4326 = {
7+
$schema: "https://proj.org/schemas/v0.7/projjson.schema.json",
8+
type: "GeographicCRS",
9+
name: "WGS 84",
10+
datum_ensemble: {
11+
name: "World Geodetic System 1984 ensemble",
12+
members: [
13+
{
14+
name: "World Geodetic System 1984 (Transit)",
15+
id: { authority: "EPSG", code: 1166 },
16+
},
17+
{
18+
name: "World Geodetic System 1984 (G730)",
19+
id: { authority: "EPSG", code: 1152 },
20+
},
21+
{
22+
name: "World Geodetic System 1984 (G873)",
23+
id: { authority: "EPSG", code: 1153 },
24+
},
25+
{
26+
name: "World Geodetic System 1984 (G1150)",
27+
id: { authority: "EPSG", code: 1154 },
28+
},
29+
{
30+
name: "World Geodetic System 1984 (G1674)",
31+
id: { authority: "EPSG", code: 1155 },
32+
},
33+
{
34+
name: "World Geodetic System 1984 (G1762)",
35+
id: { authority: "EPSG", code: 1156 },
36+
},
37+
{
38+
name: "World Geodetic System 1984 (G2139)",
39+
id: { authority: "EPSG", code: 1309 },
40+
},
41+
{
42+
name: "World Geodetic System 1984 (G2296)",
43+
id: { authority: "EPSG", code: 1383 },
44+
},
45+
],
46+
ellipsoid: {
47+
name: "WGS 84",
48+
semi_major_axis: 6378137,
49+
inverse_flattening: 298.257223563,
50+
},
51+
accuracy: "2.0",
52+
id: { authority: "EPSG", code: 6326 },
53+
},
54+
coordinate_system: {
55+
subtype: "ellipsoidal",
56+
axis: [
57+
{
58+
name: "Geodetic latitude",
59+
abbreviation: "Lat",
60+
direction: "north",
61+
unit: "degree",
62+
},
63+
{
64+
name: "Geodetic longitude",
65+
abbreviation: "Lon",
66+
direction: "east",
67+
unit: "degree",
68+
},
69+
],
70+
},
71+
scope: "Horizontal component of 3D system.",
72+
area: "World.",
73+
bbox: {
74+
south_latitude: -90,
75+
west_longitude: -180,
76+
north_latitude: 90,
77+
east_longitude: 180,
78+
},
79+
id: { authority: "EPSG", code: 4326 },
80+
};
81+
82+
describe("test TMS", () => {
83+
it("can generate TMS from EPSG CRS", async () => {
84+
const geotiff = await loadGeoTIFF(
85+
"uint8_rgb_deflate_block64_cog",
86+
"rasterio",
87+
);
88+
const crs = geotiff.crs;
89+
expect(crs).toEqual(4326);
90+
91+
const parsedCrs = wktParser(EPSG_4326);
92+
93+
const tms = generateTileMatrixSet(geotiff, parsedCrs, { id: "test-tms" });
94+
95+
expect(tms.crs).toEqual({
96+
uri: "http://www.opengis.net/def/crs/EPSG/0/4326",
97+
});
98+
expect(tms.boundingBox).toEqual({
99+
lowerLeft: [0.0, -1.28],
100+
upperRight: [1.28, 0.0],
101+
crs: { uri: "http://www.opengis.net/def/crs/EPSG/0/4326" },
102+
});
103+
expect(tms.tileMatrices).toEqual([
104+
{
105+
id: "0",
106+
scaleDenominator: 7951392.199519542,
107+
cellSize: 0.02,
108+
cornerOfOrigin: "topLeft",
109+
pointOfOrigin: [0.0, 0.0],
110+
tileWidth: 64,
111+
tileHeight: 64,
112+
matrixWidth: 1,
113+
matrixHeight: 1,
114+
},
115+
{
116+
id: "1",
117+
scaleDenominator: 3975696.099759771,
118+
cellSize: 0.01,
119+
cornerOfOrigin: "topLeft",
120+
pointOfOrigin: [0.0, 0.0],
121+
tileWidth: 64,
122+
tileHeight: 64,
123+
matrixWidth: 2,
124+
matrixHeight: 2,
125+
},
126+
]);
127+
});
128+
});

0 commit comments

Comments
 (0)