Phase 8: More maps (water-bearing sand, water level change)

app/src-tauri/src/data.rs

@@ -383,6 +383,18 @@ mod tests {
                             ..Default::default()
                         },
                     )?;
+                    // Every third one was also read three years earlier, 6 ft higher: a falling water level.
+                    if i % 3 == 0 {
+                        repo::layers::add_water_reading(
+                            tx,
+                            &b.id,
+                            &WaterReading {
+                                measured_on: format!("2023-0{}-15", 1 + i % 9),
+                                static_level: Some(44.0 + 2.0 * (i % 7) as f64),
+                                ..Default::default()
+                            },
+                        )?;
+                    }
                 }
             }
             Ok(())

app/src/model/cityMap.test.ts

@@ -1,6 +1,7 @@
 import { describe, expect, it } from "vitest";
 import { DEFAULT_MATERIALS } from "@strata/core";
-import { cityArea, cityPoint, paint, type CityGrid } from "./cityMap";
+import { THICK_RAMP, TREND_RAMP, cityArea, cityPoint, paint, ramp, thickness, type CityGrid } from "./cityMap";
+import type { Trend } from "./waterTrend";
 import { prepare, type ModelSite } from "./model";
 import { CENTRE, east, flatGrid, north, site } from "./testGround";
 
@@ -59,8 +60,8 @@ describe("cityArea", () => {
 });
 
 describe("paint", () => {
-  const grid = (codes: number[], firstWet: { top: number; level: "likely" | "rough" } | null = null): CityGrid =>
-    ({ south: 0, west: 0, north: 1, east: 1, cols: 1, rows: 1, points: [{ codes: Uint8Array.from(codes), firstWet }] });
+  const grid = (codes: number[], firstWet: { top: number; level: "likely" | "rough" } | null = null, trend: Trend | null = null): CityGrid =>
+    ({ south: 0, west: 0, north: 1, east: 1, cols: 1, rows: 1, points: [{ codes: Uint8Array.from(codes), firstWet }], trends: [trend] });
 
   it("colours the ground at a depth, paler where it is a rough guess, and leaves gaps empty", () => {
     expect([...paint(grid([1, 2, 0]), "ground", 0)]).toEqual([0xa0, 0x78, 0x5a, 200]);
@@ -78,4 +79,22 @@ describe("paint", () => {
     expect(paint(grid([5]), "thin", 0)[3]).toBe(200);
     expect(paint(grid([6]), "thin", 0)[3]).toBe(90);
   });
+
+  it("adds up the water-bearing sand down to a depth", () => {
+    // Clay, two steps of water-bearing sand (one a rough guess), clay, water-bearing sand.
+    const codes = [1, 5, 6, 1, 5];
+    expect(thickness(Uint8Array.from(codes), 20)).toEqual({ feet: 10, rough: true });
+    expect(thickness(Uint8Array.from(codes), 10)).toEqual({ feet: 5, rough: false });
+    expect(thickness(Uint8Array.from(codes), 30)).toBeNull(); // the column ends before 30 ft
+    expect(thickness(Uint8Array.from([0, 0, 0, 5]), 20)).toBeNull(); // mostly not enough data
+    expect(thickness(Uint8Array.from([1, 0, 1, 1]), 20)).toEqual({ feet: 0, rough: true }); // a gap could hide sand
+    expect([...paint(grid(codes), "thick", 10)]).toEqual([...ramp(THICK_RAMP, 5), 200]);
+    expect(paint(grid(codes), "thick", 20)[3]).toBe(90);
+  });
+
+  it("colours how fast the water level changes, pale where it is from different borewells", () => {
+    expect([...paint(grid([], null, { ftPerYear: 2, from: "same" }), "trend", 0)]).toEqual([...ramp(TREND_RAMP, 2), 200]);
+    expect(paint(grid([], null, { ftPerYear: 2, from: "area" }), "trend", 0)[3]).toBe(90);
+    expect(paint(grid([], null, null), "trend", 0)[3]).toBe(0);
+  });
 });

app/src/model/cityMap.ts

@@ -7,6 +7,7 @@ import { GROUPS } from "./groups";
 import { STEP_FT } from "./constants";
 import type { ElevationGrid } from "@strata/core";
 import { predictColumn, type Prepared } from "./model";
+import type { Trend } from "./waterTrend";
 import { GROUP_COLOUR } from "./words";
 
 export interface CityPoint {
@@ -25,9 +26,11 @@ export interface CityGrid {
   rows: number;
   /** Row by row from the north-west corner; null where the ground height is not known. */
   points: (CityPoint | null)[];
+  /** How fast the water level changes at each point, in the same order. */
+  trends: (Trend | null)[];
 }
 
-export type CityMode = "ground" | "wet" | "thin";
+export type CityMode = "ground" | "wet" | "thick" | "thin" | "trend";
 
 /** The thinnest water-bearing sand counted, as for slotted pipe. */
 const MIN_WET_FT = 10;
@@ -74,23 +77,58 @@ const rgb = (hex: string) => [1, 3, 5].map((i) => parseInt(hex.slice(i, i + 2),
 const GROUP_RGB = GROUPS.map((g) => rgb(GROUP_COLOUR[g]));
 const SURE = 200, ROUGH = 90;
 
+type Ramp = [number, number[]][];
 /** Depth to water-bearing sand: light where shallow, dark where deep. */
-export const WET_RAMP: [number, number[]][] = [[0, [207, 238, 240]], [60, [116, 186, 214]], [120, [46, 102, 179]], [200, [61, 44, 141]], [300, [30, 20, 70]]];
-export function wetColour(feet: number) {
-  const r = WET_RAMP, i = r.findIndex(([f]) => f > feet);
+export const WET_RAMP: Ramp = [[0, [207, 238, 240]], [60, [116, 186, 214]], [120, [46, 102, 179]], [200, [61, 44, 141]], [300, [30, 20, 70]]];
+/** Feet of water-bearing sand: pale where little, the sand's own colour and darker where much. */
+export const THICK_RAMP: Ramp = [[0, [247, 240, 220]], [25, [226, 196, 120]], [50, [196, 154, 60]], [100, [140, 100, 30]], [200, [80, 55, 15]]];
+/** Feet a year the water gets deeper: blue where it is rising, red where it is falling. */
+export const TREND_RAMP: Ramp = [[-5, [33, 102, 172]], [-1, [146, 197, 222]], [0, [240, 240, 240]], [1, [244, 165, 130]], [5, [178, 24, 43]]];
+
+/** The colour for a value, between the ramp's stops; past either end, the end's colour. */
+export function ramp(r: Ramp, value: number) {
+  const i = r.findIndex(([f]) => f > value);
   if (i <= 0) return i === 0 ? r[0][1] : r[r.length - 1][1];
-  const [[f0, c0], [f1, c1]] = [r[i - 1], r[i]], t = (feet - f0) / (f1 - f0);
+  const [[f0, c0], [f1, c1]] = [r[i - 1], r[i]], t = (value - f0) / (f1 - f0);
   return c0.map((v, k) => Math.round(v + (c1[k] - v) * t));
 }
+
+const WET = GROUPS.indexOf("waterBearing");
+
+/**
+ * Feet of water-bearing sand proposed from the surface down to `depthFt`, and whether any of it is
+ * a rough guess. Null where the column ends before that depth or more than half of it has no proposal.
+ */
+export function thickness(codes: Uint8Array, depthFt: number) {
+  const steps = Math.floor(depthFt / STEP_FT);
+  if (steps === 0 || codes.length < steps) return null;
+  let wet = 0, rough = false, unknown = 0;
+  for (let i = 0; i < steps; i++) {
+    const code = codes[i];
+    if (!code) unknown++;
+    else if ((code - 1) >> 1 === WET) {
+      wet++;
+      rough ||= ((code - 1) & 1) === 1;
+    }
+  }
+  // A gap in the proposal could hide sand, so the total is then only a rough guess.
+  return unknown * 2 > steps ? null : { feet: wet * STEP_FT, rough: rough || unknown > 0 };
+}
 const THIN_RGB = [13, 104, 131];
 
-/** The map as RGBA pixels, one per point, rows from the north. `depthFt` is used by "ground" and "thin". */
+/** The map as RGBA pixels, one per point, rows from the north. `depthFt` is used by "ground", "thick" and "thin". */
 export function paint(g: CityGrid, mode: CityMode, depthFt: number): Uint8ClampedArray<ArrayBuffer> {
   const out = new Uint8ClampedArray(g.points.length * 4), step = Math.floor(depthFt / STEP_FT);
   g.points.forEach((p, i) => {
     let colour: number[] | null = null, alpha = 0;
-    if (p && mode === "wet") {
-      if (p.firstWet) [colour, alpha] = [wetColour(p.firstWet.top), p.firstWet.level === "likely" ? SURE : ROUGH];
+    const trend = g.trends[i];
+    if (mode === "trend") {
+      if (trend) [colour, alpha] = [ramp(TREND_RAMP, trend.ftPerYear), trend.from === "same" ? SURE : ROUGH];
+    } else if (p && mode === "wet") {
+      if (p.firstWet) [colour, alpha] = [ramp(WET_RAMP, p.firstWet.top), p.firstWet.level === "likely" ? SURE : ROUGH];
+    } else if (p && mode === "thick") {
+      const t = thickness(p.codes, depthFt);
+      if (t) [colour, alpha] = [ramp(THICK_RAMP, t.feet), t.rough ? ROUGH : SURE];
     } else if (p) {
       const code = p.codes[step] ?? 0;
       if (code) [colour, alpha] = [mode === "ground" ? GROUP_RGB[(code - 1) >> 1] : THIN_RGB, (code - 1) & 1 ? ROUGH : SURE];

app/src/model/waterTrend.test.ts

@@ -0,0 +1,78 @@
+import { describe, expect, it } from "vitest";
+import type { Borewell, WaterReading } from "@strata/core";
+import { lineThrough, readingsOf, sameBorewellTrends, trendAt } from "./waterTrend";
+import { CENTRE, east, north } from "./testGround";
+
+const bw = (id: string, km: number, level: number | null = null, on: string | null = null) =>
+  ({ id, borewellId: id.toUpperCase(), latitude: north(0), longitude: east(km), waterLevel: level, waterLevelOn: on }) as Borewell;
+const reading = (measuredOn: string, staticLevel: number | null) => ({ measuredOn, staticLevel }) as WaterReading;
+const site = (b: Borewell, readings: WaterReading[] = []) => ({ borewell: b, waterReadings: readings });
+
+describe("lineThrough", () => {
+  it("gives the change per year of a straight line through the readings", () => {
+    expect(lineThrough([{ year: 2020, feet: 40 }, { year: 2022, feet: 46 }, { year: 2024, feet: 52 }])).toBeCloseTo(3, 5);
+  });
+});
+
+describe("readingsOf", () => {
+  it("takes dated readings, or a borewell's own dated level when it has none", () => {
+    const r = readingsOf([
+      site(bw("a", 0, 50, "2025-06-01"), [reading("2023-01-01", 40), reading("2025-01-01", 44), reading("2024-01-01", null)]),
+      site(bw("b", 1, 60, "2024-07-01")),
+      site(bw("c", 2, 70, null)),
+    ]);
+    expect(r.map((x) => [x.id, x.feet])).toEqual([["a", 40], ["a", 44], ["a", 50], ["b", 60]]);
+    expect(r[0].year).toBeCloseTo(2023, 2);
+  });
+
+  it("keeps a borewell's own level when its readings only have the level with the pump running", () => {
+    const r = readingsOf([site(bw("a", 0, 50, "2023-06-01"), [{ measuredOn: "2026-01-01", staticLevel: null, dynamicLevel: 80 } as WaterReading])]);
+    expect(r.map((x) => x.feet)).toEqual([50]);
+  });
+
+  it("does not count the newest reading twice", () => {
+    const r = readingsOf([site(bw("a", 0, 44, "2025-01-01"), [reading("2023-01-01", 40), reading("2025-01-01", 44)])]);
+    expect(r).toHaveLength(2);
+  });
+
+  it("reads dates with spaces around them", () => {
+    expect(readingsOf([site(bw("a", 0, 50, " 2024-01-01 "))])).toHaveLength(1);
+  });
+});
+
+describe("trendAt", () => {
+  it("prefers borewells measured more than once, nearer ones counting more", () => {
+    const readings = readingsOf([
+      site(bw("near", 0.2), [reading("2022-01-01", 40), reading("2024-01-01", 44)]), // +2 a year
+      site(bw("far", 1.5), [reading("2022-01-01", 40), reading("2024-01-01", 52)]), // +6 a year
+      ...Array.from({ length: 6 }, (_, i) => site(bw(`o${i}`, 0.5, 30, `${2020 + i}-01-01`))),
+    ]);
+    const t = trendAt(CENTRE[0], CENTRE[1], readings, sameBorewellTrends(readings))!;
+    expect(t.from).toBe("same");
+    expect(t.ftPerYear).toBeGreaterThan(2);
+    expect(t.ftPerYear).toBeLessThan(2.5);
+  });
+
+  it("ignores repeat readings less than a year apart", () => {
+    const readings = readingsOf([site(bw("a", 0.2), [reading("2024-01-01", 40), reading("2024-06-01", 48)])]);
+    expect(sameBorewellTrends(readings)).toEqual([]);
+  });
+
+  it("otherwise uses different borewells over the years, if there are enough", () => {
+    const readings = readingsOf(Array.from({ length: 15 }, (_, i) => site(bw(`o${i}`, 0.1 * i, 40 + (i % 5), `${2021 + (i % 5)}-01-01`))));
+    expect(trendAt(CENTRE[0], CENTRE[1], readings, [])).toMatchObject({ from: "area", ftPerYear: expect.closeTo(1, 3) });
+  });
+
+  it("says nothing with fewer than 15 readings or under 2 years", () => {
+    // Different wells differ by far more than a few years' change, so a handful of them says nothing.
+    const few = readingsOf(Array.from({ length: 14 }, (_, i) => site(bw(`o${i}`, 0.1 * i, 40 + i, `${2021 + (i % 5)}-01-01`))));
+    expect(trendAt(CENTRE[0], CENTRE[1], few, [])).toBeNull();
+    const oneYear = readingsOf(Array.from({ length: 16 }, (_, i) => site(bw(`o${i}`, 0.1 * i, 40 + i, `2024-0${1 + (i % 9)}-01`))));
+    expect(trendAt(CENTRE[0], CENTRE[1], oneYear, [])).toBeNull();
+  });
+
+  it("looks only within 2 km", () => {
+    const readings = readingsOf(Array.from({ length: 6 }, (_, i) => site(bw(`o${i}`, 3 + i * 0.1, 40 + i, `${2020 + i}-01-01`))));
+    expect(trendAt(CENTRE[0], CENTRE[1], readings, sameBorewellTrends(readings))).toBeNull();
+  });
+});

app/src/model/waterTrend.ts

@@ -0,0 +1,98 @@
+/**
+ * How fast the water level is changing: feet per year it gets deeper (positive) or shallower.
+ * Measured at the same borewells where they were read more than once, a year or more apart;
+ * otherwise, roughly, from different borewells nearby measured in different years.
+ */
+import type { Borewell, WaterReading } from "@strata/core";
+
+export interface Reading {
+  /** The borewell's record id. */
+  id: string;
+  latitude: number;
+  longitude: number;
+  /** When, as a year with its fraction (2024.5 is about the start of July 2024). */
+  year: number;
+  /** Feet below ground. */
+  feet: number;
+}
+
+export interface Trend {
+  ftPerYear: number;
+  from: "same" | "area";
+}
+
+/** Within this distance of a square its readings count. */
+const NEAR_KM = 2;
+/** Distances closer than this count as this, so one borewell on the square does not count without limit. */
+const CLOSEST_KM = 0.1;
+/**
+ * The area trend compares different wells, which differ by far more than a few years' change; with
+ * few readings it would mostly show that difference as a strong trend. So it needs many.
+ */
+const AREA = { readings: 15, years: 2 };
+
+const yearOf = (iso: string) => {
+  const d = new Date(`${iso.trim().slice(0, 10)}T00:00:00Z`);
+  const y = d.getUTCFullYear();
+  return Number.isNaN(y) ? null : y + (d.getTime() - Date.UTC(y, 0, 1)) / (365.25 * 86_400_000);
+};
+
+/**
+ * Every dated water level of located borewells: their readings with the level at rest, and their own
+ * level too unless a reading of that day already gives it (it is usually the newest reading's).
+ */
+export function readingsOf(sites: { borewell: Borewell; waterReadings: WaterReading[] }[]): Reading[] {
+  const out: Reading[] = [];
+  for (const { borewell: b, waterReadings } of sites) {
+    if (b.latitude == null || b.longitude == null) continue;
+    const dated = waterReadings.filter((r) => r.staticLevel != null).map((r) => [r.measuredOn, r.staticLevel!] as const);
+    const day = (on: string) => on.trim().slice(0, 10);
+    if (b.waterLevelOn && b.waterLevel != null && !dated.some(([on]) => day(on) === day(b.waterLevelOn!))) dated.push([b.waterLevelOn, b.waterLevel]);
+    for (const [on, feet] of dated) {
+      const year = yearOf(on);
+      if (year != null) out.push({ id: b.id, latitude: b.latitude, longitude: b.longitude, year, feet });
+    }
+  }
+  return out;
+}
+
+/** Feet per year of the straight line through the readings (least squares); null without a spread of dates. */
+export function lineThrough(points: { year: number; feet: number }[]): number | null {
+  const n = points.length, my = points.reduce((s, p) => s + p.year, 0) / n, mf = points.reduce((s, p) => s + p.feet, 0) / n;
+  let sxy = 0, sxx = 0;
+  for (const p of points) {
+    sxy += (p.year - my) * (p.feet - mf);
+    sxx += (p.year - my) ** 2;
+  }
+  return sxx > 0 ? sxy / sxx : null;
+}
+
+const span = (rs: Reading[]) => Math.max(...rs.map((r) => r.year)) - Math.min(...rs.map((r) => r.year));
+
+/** The change per year at each borewell read at least twice, a year or more apart. */
+export function sameBorewellTrends(readings: Reading[]): { latitude: number; longitude: number; ftPerYear: number }[] {
+  const byBorewell = new Map<string, Reading[]>();
+  for (const r of readings) byBorewell.set(r.id, [...(byBorewell.get(r.id) ?? []), r]);
+  return [...byBorewell.values()]
+    .filter((rs) => rs.length >= 2 && span(rs) >= 1)
+    .map((rs) => ({ latitude: rs[0].latitude, longitude: rs[0].longitude, ftPerYear: lineThrough(rs)! }));
+}
+
+export function trendAt(latitude: number, longitude: number, readings: Reading[], same: ReturnType<typeof sameBorewellTrends>): Trend | null {
+  const kx = 111.32 * Math.cos((latitude * Math.PI) / 180);
+  const km = (p: { latitude: number; longitude: number }) => Math.hypot((p.longitude - longitude) * kx, (p.latitude - latitude) * 110.57);
+
+  let total = 0, weighted = 0;
+  for (const t of same) {
+    const d = km(t);
+    if (d > NEAR_KM) continue;
+    const w = 1 / Math.max(d, CLOSEST_KM) ** 2;
+    total += w;
+    weighted += w * t.ftPerYear;
+  }
+  if (total > 0) return { ftPerYear: weighted / total, from: "same" };
+
+  const near = readings.filter((r) => km(r) <= NEAR_KM);
+  if (near.length < AREA.readings || span(near) < AREA.years) return null;
+  return { ftPerYear: lineThrough(near)!, from: "area" };
+}

app/src/text/index.ts

@@ -115,20 +115,31 @@ export const text = {
 
   city: {
     modeLabel: "What the map shows",
-    modes: { ground: "Ground at a depth", wet: "Depth to water-bearing sand", thin: "Where data is thin" },
+    modes: {
+      ground: "Ground at a depth", wet: "Depth to water-bearing sand", thick: "Water-bearing sand down to a depth",
+      thin: "Where data is thin", trend: "How fast the water level is changing",
+    },
     about: {
       ground: "What the ground is most likely made of at the chosen depth below ground.",
       wet: "How deep the first water-bearing sand at least 10 ft thick is likely to be.",
+      thick: "How many feet of water-bearing sand are likely between the surface and the chosen depth.",
       thin: "How sure the proposals are at the chosen depth, on their own, so the gaps in the records stand out.",
+      trend: "How many feet a year the water level is getting deeper (falling) or shallower (rising), from water levels measured within 2 km.",
     },
     depth: (feet: number) => `Depth: ${feet} ft below ground`,
     working: (percent: number) => `Working out the map… ${percent}%`,
-    wetKey: "Feet below ground",
+    rampKey: { wet: "Feet below ground", thick: "Feet of water-bearing sand", trend: "Feet a year deeper (+) or shallower (−): red where the water is falling, blue where it is rising" },
+    strongTrend: "Strong: measured more than once at the same borewells",
+    palerTrend: "Paler: from different borewells measured in different years nearby (a rough guess)",
+    repeat: (n: number) =>
+      `${n === 1 ? "1 borewell has" : `${n} borewells have`} water levels measured a year or more apart. The map firms up as more repeat readings are added in StrataField.`,
     paler: "Paler: a rough guess",
     clear: {
       ground: "Clear: not enough data at this depth",
       wet: "Clear: no such sand proposed, or not enough data",
+      thick: "Clear: not enough data down to this depth",
       thin: "Clear: not enough data at this depth",
+      trend: "Clear: too few dated water levels within 2 km (at least 15, over 2 years or more)",
     },
     squares: (metres: number) => `Each square is about ${metres} m across.`,
     empty: "There are no borewells with a location and layers to map yet.",