diff --git a/apps/class-solid/src/components/Analysis.tsx b/apps/class-solid/src/components/Analysis.tsx index 3b0b5269..8ac6f291 100644 --- a/apps/class-solid/src/components/Analysis.tsx +++ b/apps/class-solid/src/components/Analysis.tsx @@ -1,8 +1,9 @@ import { For, Match, Show, Switch, createMemo, createUniqueId } from "solid-js"; -import { getVerticalProfiles } from "~/lib/profiles"; +import { getThermodynamicProfiles, getVerticalProfiles } from "~/lib/profiles"; import { type Analysis, deleteAnalysis, experiments } from "~/lib/store"; -import LinePlot from "./LinePlot"; import { MdiCog, MdiContentCopy, MdiDelete, MdiDownload } from "./icons"; +import LinePlot from "./plots/LinePlot"; +import { SkewTPlot } from "./plots/skewTlogP"; import { Button } from "./ui/button"; import { Card, CardContent, CardHeader, CardTitle } from "./ui/card"; @@ -37,19 +38,25 @@ export function TimeSeriesPlot() { .map((perm, j) => { return { label: `${e.name}/${perm.name}`, - y: perm.output?.h ?? [], - x: perm.output?.t ?? [], color: colors[(j + 1) % 10], linestyle: linestyles[i % 5], + data: + perm.output?.t.map((tVal, i) => ({ + x: tVal, + y: perm.output?.h[i] || Number.NaN, + })) || [], }; }); return [ { - y: experimentOutput?.h ?? [], - x: experimentOutput?.t ?? [], label: e.name, color: colors[0], linestyle: linestyles[i], + data: + experimentOutput?.t.map((tVal, i) => ({ + x: tVal, + y: experimentOutput?.h[i] || Number.NaN, + })) || [], }, ...permutationRuns, ]; @@ -79,7 +86,7 @@ export function VerticalProfilePlot() { color: colors[(j + 1) % 10], linestyle: linestyles[i % 5], label: `${e.name}/${p.name}`, - ...getVerticalProfiles(p.output, p.config, variable, time), + data: getVerticalProfiles(p.output, p.config, variable, time), }; }); @@ -88,7 +95,7 @@ export function VerticalProfilePlot() { label: e.name, color: colors[0], linestyle: linestyles[i], - ...getVerticalProfiles( + data: getVerticalProfiles( e.reference.output ?? { t: [], h: [], @@ -113,6 +120,39 @@ export function VerticalProfilePlot() { ); } +export function ThermodynamicPlot() { + const time = -1; + const skewTData = createMemo(() => { + return experiments.flatMap((e, i) => { + const permutations = e.permutations.map((p, j) => { + // TODO get additional config info from reference + // permutations probably usually don't have gammaq/gammatetha set? + return { + color: colors[(j + 1) % 10], + linestyle: linestyles[i % 5], + label: `${e.name}/${p.name}`, + data: getThermodynamicProfiles(p.output, p.config, time), + }; + }); + + return [ + { + label: e.name, + color: colors[0], + linestyle: linestyles[i], + data: getThermodynamicProfiles( + e.reference.output, + e.reference.config, + time, + ), + }, + ...permutations, + ]; + }); + }); + return ; +} + /** Simply show the final height for each experiment that has output */ function FinalHeights() { return ( @@ -154,7 +194,7 @@ function FinalHeights() { export function AnalysisCard(analysis: Analysis) { const id = createUniqueId(); return ( - + {/* TODO: make name & description editable */} {analysis.name} @@ -183,6 +223,7 @@ export function AnalysisCard(analysis: Analysis) { Unknown analysis type

}> + {/* @ts-ignore: kept for developers, but not included in production */} @@ -192,6 +233,9 @@ export function AnalysisCard(analysis: Analysis) { + + +
diff --git a/apps/class-solid/src/components/LinePlot.tsx b/apps/class-solid/src/components/LinePlot.tsx deleted file mode 100644 index 8916bad6..00000000 --- a/apps/class-solid/src/components/LinePlot.tsx +++ /dev/null @@ -1,126 +0,0 @@ -import * as d3 from "d3"; -import { For } from "solid-js"; -import { AxisBottom, AxisLeft } from "./Axes"; - -export interface ChartData { - label: string; - color: string; - linestyle: string; - x: number[]; - y: number[]; -} - -/** - * Calculate a "nice" step size by rounding up to the nearest power of 10 - * Snap the min and max to the nearest multiple of step - */ -function getNiceAxisLimits(data: number[]): [number, number] { - const max = Math.max(...data); - const min = Math.min(...data); - const range = max - min; - const step = 10 ** Math.floor(Math.log10(range)); - - const niceMin = Math.floor(min / step) * step; - const niceMax = Math.ceil(max / step) * step; - - return [niceMin, niceMax]; -} - -/** - * Zip x and y such that we can iterate over pairs of [x, y] in d3.line - */ -function zipXY(data: ChartData): [number, number][] { - const length = data.x.length; - return Array.from({ length }, (_, i) => [data.x[i], data.y[i]]); -} - -export default function LinePlot({ - data, - xlabel, - ylabel, -}: { data: () => ChartData[]; xlabel?: string; ylabel?: string }) { - // TODO: Make responsive - const width = 450; - const height = 400; - const [marginTop, marginRight, marginBottom, marginLeft] = [25, 50, 50, 50]; - - const xLim = () => getNiceAxisLimits(data().flatMap((d) => d.x)); - const yLim = () => getNiceAxisLimits(data().flatMap((d) => d.y)); - const scaleX = () => - d3.scaleLinear(xLim(), [marginLeft, width - marginRight]); - const scaleY = () => - d3.scaleLinear(yLim(), [height - marginBottom, marginTop]); - - // const l = d3.line((d, i) => scaleX(x[i]), scaleY); - const l = d3.line( - (d) => scaleX()(d[0]), - (d) => scaleY()(d[1]), - ); - - return ( -
- {/* Legend */} -
- - {(d) => ( - <> - - - legend - - -

{d.label}

-
- - )} -
-
- - {/* Plot */} - - Vertical profile plot - {/* Axes */} - - - - {/* Line */} - - {(d) => ( - - {d.label} - - )} - - -
- ); -} diff --git a/apps/class-solid/src/components/Axes.tsx b/apps/class-solid/src/components/plots/Axes.tsx similarity index 52% rename from apps/class-solid/src/components/Axes.tsx rename to apps/class-solid/src/components/plots/Axes.tsx index 6b3f49eb..f0f8b480 100644 --- a/apps/class-solid/src/components/Axes.tsx +++ b/apps/class-solid/src/components/plots/Axes.tsx @@ -9,23 +9,26 @@ interface AxisProps { transform?: string; tickCount?: number; label?: string; + tickValues?: number[]; + tickFormat?: (n: number | { valueOf(): number }) => string; + decreasing?: boolean; } -export const AxisBottom = (props: AxisProps) => { - const ticks = () => { - const domain = props.scale.domain(); - const tickCount = props.tickCount || 5; +const ticks = (props: AxisProps) => { + const domain = props.scale.domain(); + const generateTicks = (domain = [0, 1], tickCount = 5) => { const step = (domain[1] - domain[0]) / (tickCount - 1); - - return [...Array(tickCount).keys()].map((i) => { - const value = domain[0] + i * step; - return { - value, - position: props.scale(value), - }; - }); + return [...Array(10).keys()].map((i) => domain[0] + i * step); }; + const values = props.tickValues + ? props.tickValues.filter((x) => x >= domain[0] && x <= domain[1]) + : generateTicks(domain, props.tickCount); + return values.map((value) => ({ value, position: props.scale(value) })); +}; + +export const AxisBottom = (props: AxisProps) => { + const format = props.tickFormat ? props.tickFormat : d3.format(".3g"); return ( { y2="0" stroke="currentColor" /> - + {(tick) => ( - {d3.format(".3g")(tick.value)} + {format(tick.value)} )} @@ -53,21 +56,8 @@ export const AxisBottom = (props: AxisProps) => { }; export const AxisLeft = (props: AxisProps) => { - const ticks = () => { - const domain = props.scale.domain(); - const tickCount = props.tickCount || 5; - const step = (domain[1] - domain[0]) / (tickCount - 1); - - return [...Array(tickCount).keys()].map((i) => { - const value = domain[0] + i * step; - return { - value, - position: props.scale(value), - }; - }); - }; - - const labelpos = props.scale.range().reduce((a, b) => a + b) / 2; + const format = props.tickFormat ? props.tickFormat : d3.format(".0f"); + const yAnchor = props.decreasing ? 0 : 1; return ( { y2={props.scale.range()[1]} stroke="currentColor" /> - + {(tick) => ( - {tick.value.toFixed()} + {format(tick.value)} )} {props.label} ); }; + +/** + * Calculate a "nice" step size by rounding up to the nearest power of 10 + * Snap the min and max to the nearest multiple of step + */ +export function getNiceAxisLimits(data: number[]): [number, number] { + const max = Math.max(...data); + const min = Math.min(...data); + const range = max - min; + const step = 10 ** Math.floor(Math.log10(range)); + + const niceMin = Math.floor(min / step) * step; + const niceMax = Math.ceil(max / step) * step; + + return [niceMin, niceMax]; +} diff --git a/apps/class-solid/src/components/plots/Base.tsx b/apps/class-solid/src/components/plots/Base.tsx new file mode 100644 index 00000000..d0aef4e1 --- /dev/null +++ b/apps/class-solid/src/components/plots/Base.tsx @@ -0,0 +1,9 @@ +export interface ChartData { + label: string; + color: string; + linestyle: string; + data: T[]; +} + +// TODO: would be nice to create a chartContainer/context that manages logic like +// width/height/margins etc. that should be consistent across different plots. diff --git a/apps/class-solid/src/components/plots/Legend.tsx b/apps/class-solid/src/components/plots/Legend.tsx new file mode 100644 index 00000000..c3be8ad2 --- /dev/null +++ b/apps/class-solid/src/components/plots/Legend.tsx @@ -0,0 +1,45 @@ +import { For } from "solid-js"; +import { cn } from "~/lib/utils"; +import type { ChartData } from "./Base"; + +export interface LegendProps { + entries: () => ChartData[]; + width: string; +} + +export function Legend(props: LegendProps) { + return ( + // {/* Legend */} +
+ + {(d) => ( + <> + + + legend + + +

{d.label}

+
+ + )} +
+
+ ); +} diff --git a/apps/class-solid/src/components/plots/LinePlot.tsx b/apps/class-solid/src/components/plots/LinePlot.tsx new file mode 100644 index 00000000..d18a2438 --- /dev/null +++ b/apps/class-solid/src/components/plots/LinePlot.tsx @@ -0,0 +1,78 @@ +import * as d3 from "d3"; +import { For } from "solid-js"; +import { AxisBottom, AxisLeft, getNiceAxisLimits } from "./Axes"; +import type { ChartData } from "./Base"; +import { Legend } from "./Legend"; + +export interface Point { + x: number; + y: number; +} + +export default function LinePlot({ + data, + xlabel, + ylabel, +}: { data: () => ChartData[]; xlabel?: string; ylabel?: string }) { + // TODO: Make responsive + // const margin = [30, 40, 20, 45]; // reference from skew-T + const [marginTop, marginRight, marginBottom, marginLeft] = [20, 20, 35, 55]; + const width = 500; + const height = 500; + const w = 500 - marginRight - marginLeft; + const h = 500 - marginTop - marginBottom; + + const xLim = () => + getNiceAxisLimits(data().flatMap((d) => d.data.flatMap((d) => d.x))); + const yLim = () => + getNiceAxisLimits(data().flatMap((d) => d.data.flatMap((d) => d.y))); + const scaleX = () => d3.scaleLinear(xLim(), [0, w]); + const scaleY = () => d3.scaleLinear(yLim(), [h, 0]); + + const l = d3.line( + (d) => scaleX()(d.x), + (d) => scaleY()(d.y), + ); + + return ( +
+ + {/* Plot */} + + + Vertical profile plot + {/* Axes */} + + + + {/* Line */} + + {(d) => ( + + {d.label} + + )} + + + +
+ ); +} diff --git a/apps/class-solid/src/components/plots/skewTlogP.tsx b/apps/class-solid/src/components/plots/skewTlogP.tsx new file mode 100644 index 00000000..e9de0430 --- /dev/null +++ b/apps/class-solid/src/components/plots/skewTlogP.tsx @@ -0,0 +1,198 @@ +// Code modified from https://github.com/rsobash/d3-skewt/ (MIT license) +import * as d3 from "d3"; +import { For, createSignal } from "solid-js"; +import { AxisBottom, AxisLeft } from "./Axes"; +import type { ChartData } from "./Base"; +import { Legend } from "./Legend"; + +interface SoundingRecord { + p: number; + T: number; + Td: number; +} + +const deg2rad = Math.PI / 180; +const tan = Math.tan(55 * deg2rad); +const basep = 1050; +const topPressure = 100; + +function SkewTBackGround({ + w, + h, + x, + y, +}: { + w: number; + h: number; + x: d3.ScaleLinear; + y: d3.ScaleLogarithmic; +}) { + const pressureLines = [1000, 850, 700, 500, 300, 200, 100]; + const temperatureLines = d3.range(-100, 45, 10); + + // Dry adiabats (lines of constant potential temperature): array of lines of [p, T] + const pressurePoints = d3.range(topPressure, basep + 1, 10); + const temperaturePoints = d3.range(-30, 240, 20); + const dryAdiabats: [number, number][][] = temperaturePoints.map( + (temperature) => pressurePoints.map((pressure) => [pressure, temperature]), + ); + + const dryline = d3 + .line() + .x( + (d) => + x((273.15 + d[1]) / (1000 / d[0]) ** 0.286 - 273.15) + + (y(basep) - y(d[0])) / tan, + ) + .y((d) => y(d[0])); + return ( + + + + + {/* Add grid */} + {/* Temperature lines */} + + {(tline) => ( + + )} + + {/* Pressure lines */} + + {(pline) => ( + + )} + + {/* Dry Adiabats */} + + {(d) => ( + + )} + + {/* Line along right edge of plot */} + + + + + + ); +} + +// Note: using temperatures in Kelvin as that's easiest to get from CLASS, but +// perhaps not the most interoperable with other sounding data sources. +export function SkewTPlot({ + data, +}: { data: () => ChartData[] }) { + const [hovered, setHovered] = createSignal(null); + const width = 500; + const height = 500; + const [marginTop, marginRight, marginBottom, marginLeft] = [20, 20, 35, 55]; + const w = 500 - marginRight - marginLeft; + const h = 500 - marginTop - marginBottom; + + // Scales and axes. Note the inverted domain for the y-scale: bigger is up! + const x = d3.scaleLinear().range([0, w]).domain([-45, 50]); + const y = d3.scaleLog().range([0, h]).domain([topPressure, basep]); + + const temperatureLine = d3 + .line() + .x((d) => x(d.T - 273.15) + (y(basep) - y(d.p)) / tan) + .y((d) => y(d.p)); + + const dewpointLine = d3 + .line() + .x((d) => x(d.Td - 273.15) + (y(basep) - y(d.p)) / tan) + .y((d) => y(d.p)); + + // // bisector function for tooltips + // const bisectTemp = d3.bisector((d) => d.press).left; + + return ( +
+ + {/* Create svg container for sounding */} + + Thermodynamic diagram + + + + {(d, index) => ( + setHovered(index())} + onMouseLeave={() => setHovered(null)} + > + {d.label} + + + + )} + + + +
+ ); +} diff --git a/apps/class-solid/src/lib/profiles.ts b/apps/class-solid/src/lib/profiles.ts index d21cd3a6..6519c56c 100644 --- a/apps/class-solid/src/lib/profiles.ts +++ b/apps/class-solid/src/lib/profiles.ts @@ -1,5 +1,6 @@ import type { ClassOutput } from "@classmodel/class/runner"; import type { PartialConfig } from "@classmodel/class/validate"; +import type { Point } from "~/components/plots/LinePlot"; // Get vertical profiles for a single class run export function getVerticalProfiles( @@ -7,39 +8,152 @@ export function getVerticalProfiles( config: PartialConfig, variable: "theta" | "q" = "theta", t = -1, -): { x: number[]; y: number[] } { +): Point[] { // Guard against undefined output if (output === undefined) { - return { y: [], x: [] }; + return []; } // Extract height profile const height = output.h.slice(t)[0]; const dh = 400; // how much free troposphere to display? - const h = [0, height, height, height + dh]; + const hProfile = [0, height, height, height + dh]; if (variable === "theta") { // Extract potential temperature profile - const thetas = output.theta.slice(t)[0]; + const theta = output.theta.slice(t)[0]; const dtheta = output.dtheta.slice(t)[0]; - const gammatheta = config.mixedLayer?.gammatheta ?? 0; - const theta = [ - thetas, - thetas, - thetas + dtheta, - thetas + dtheta + dh * gammatheta, + // TODO: make sure config contains gammatheta + const gammatheta = config.mixedLayer?.gammatheta ?? 0.006; + const thetaProfile = [ + theta, + theta, + theta + dtheta, + theta + dtheta + dh * gammatheta, ]; - return { y: h, x: theta }; + return hProfile.map((h, i) => ({ x: thetaProfile[i], y: h })); } - if (variable === "q") { // Extract humidity profile - const qs = output.q.slice(t)[0]; + const q = output.q.slice(t)[0]; const dq = output.dq.slice(t)[0]; + // TODO: make sure config contains gammaq const gammaq = config.mixedLayer?.gammaq ?? 0; - const q = [qs, qs, qs + dq, qs + dq + dh * gammaq]; - return { y: h, x: q }; + const qProfile = [q, q, q + dq, q + dq + dh * gammaq]; + return hProfile.map((h, i) => ({ x: qProfile[i], y: h })); + } + return []; +} + +/** + * https://en.wikipedia.org/wiki/Dew_point#Calculating_the_dew_point + */ +const dewpoint = (q: number, p: number) => { + // Empirical fit parameters (Sonntag1990, see wikipedia entry for more options) + const A = 6.112; + const B = 17.62; + const C = 243.12; + + const w = q / (1 - q); // mixing ratio + const e = (w * p) / (w + 0.622); // Actual vapour pressure; Wallace and Hobbs 3.59 + + const Td = (C * Math.log(e / A)) / (B - Math.log(e / A)); + + return Td + 273.15; +}; + +/** + * Calculate temperature from potential temperature + * https://unidata.github.io/MetPy/latest/api/generated/metpy.calc.temperature_from_potential_temperature.html# + */ +const thetaToT = (theta: number, p: number, p0 = 1000) => { + const R = 287; // specific gas constant for dry air + const cp = 1004; // specific heat of dry air at constant pressure + const T = theta * (p / p0) ** (R / cp); + return T; +}; + +/** + * Calculate pressure difference over layer using hypsometric equation + * Wallace & Hobbs eq (3.23) + */ +const pressureDiff = (T: number, q: number, p: number, dz: number) => { + const Rd = 287; // specific gas constant for dry air + const g = 9.81; // gravity + const w = q / (1 - q); // mixing ratio + const Tv = T * (1 + 0.61 * w); // virtual temperature + const dp = (dz * p * g) / (Rd * Tv); + return -dp; +}; + +/** + * Calculate thickness of layer using hypsometric equation + * Wallace & Hobbs eq (3.23) + */ +const thickness = (T: number, q: number, p: number, dp: number) => { + const Rd = 287; // specific gas constant for dry air + const g = 9.81; // gravity + const w = q / (1 - q); // mixing ratio + const Tv = T * (1 + 0.61 * w); // virtual temperature + const dz = (dp * Rd * Tv) / (p * g); + return dz; +}; + +export function getThermodynamicProfiles( + output: ClassOutput | undefined, + config: PartialConfig, + t = -1, +) { + // Guard against undefined output + if (output === undefined) { + return []; + } + + let theta = output.theta.slice(t)[0]; + let q = output.q.slice(t)[0]; + const dtheta = output.dtheta.slice(t)[0]; + const dq = output.dq.slice(t)[0]; + const h = output.h.slice(t)[0]; + // TODO: ensure config contains gammatheta and gammaq + const gammaTheta = config.mixedLayer?.gammatheta ?? 0.006; + const gammaq = config.mixedLayer?.gammaq ?? 0; + + const nz = 25; + let dz = h / nz; + const zArrayMixedLayer = [...Array(nz).keys()].map((i) => i * dz); + + let p = 1000; // hPa?? + let T = theta; + let Td = dewpoint(q, p); + const soundingData: { p: number; T: number; Td: number }[] = [{ p, T, Td }]; + + // Mixed layer + for (const z of zArrayMixedLayer) { + p += pressureDiff(T, q, p, dz); + T = thetaToT(theta, p); + Td = dewpoint(q, p); + + soundingData.push({ p, T, Td }); + } + + // Inversion + theta += dtheta; + q += dq; + T = thetaToT(theta, p); + Td = dewpoint(q, p); + soundingData.push({ p, T, Td }); + + // Free troposphere + dz = 200; + while (p > 100) { + theta += dz * gammaTheta; + q += dz * gammaq; + p += pressureDiff(T, q, p, dz); + T = thetaToT(theta, p); + Td = dewpoint(q, p); + + soundingData.push({ p, T, Td }); } - return { y: [], x: [] }; + return soundingData; } diff --git a/apps/class-solid/src/lib/store.ts b/apps/class-solid/src/lib/store.ts index a1bb6210..94709c11 100644 --- a/apps/class-solid/src/lib/store.ts +++ b/apps/class-solid/src/lib/store.ts @@ -260,12 +260,13 @@ export async function loadStateFromString(rawState: string): Promise { await Promise.all(loadedExperiments.map((_, i) => runExperiment(i))); } -const analysisNames = { +export const analysisNames = { profiles: "Vertical profiles", timeseries: "Timeseries", - finalheight: "Final height", + skewT: "Thermodynamic diagram", + // finalheight: "Final height", // keep for development but not in production } as const; -type AnalysisType = keyof typeof analysisNames; +export type AnalysisType = keyof typeof analysisNames; export interface Analysis { name: string; diff --git a/apps/class-solid/src/routes/index.tsx b/apps/class-solid/src/routes/index.tsx index 24024ddd..39089655 100644 --- a/apps/class-solid/src/routes/index.tsx +++ b/apps/class-solid/src/routes/index.tsx @@ -16,7 +16,12 @@ import { Flex } from "~/components/ui/flex"; import { Toaster } from "~/components/ui/toast"; import { onPageLoad } from "~/lib/state"; -import { addAnalysis, experiments } from "~/lib/store"; +import { + type AnalysisType, + addAnalysis, + analysisNames, + experiments, +} from "~/lib/store"; import { analyses } from "~/lib/store"; export default function Home() { @@ -60,18 +65,15 @@ export default function Home() { Add analysis - addAnalysis("finalheight")}> - Final height - - addAnalysis("timeseries")}> - Timeseries - - addAnalysis("profiles")}> - Vertical profile - - - Skew-T diagram (not implemented) - + + {([key, value]) => ( + addAnalysis(key)}> + {value} + + )} + diff --git a/apps/class-solid/tests/share.spec.ts b/apps/class-solid/tests/share.spec.ts index eb11f88c..b5ac985c 100644 --- a/apps/class-solid/tests/share.spec.ts +++ b/apps/class-solid/tests/share.spec.ts @@ -31,18 +31,19 @@ test("Create share link from an experiment", async ({ page }) => { const config1 = await parseDownload(downloadPromise1); expect(config1.reference.initialState?.h_0).toEqual(800); - // Check that shared experiment has been run by - // checking height in final height analysis - const finalHeightAnalysis = sharedPage.getByRole("article", { - name: "Final height", - }); - const finalHeightOfExperiment = finalHeightAnalysis.getByRole("listitem", { - name: "My experiment 1", - exact: true, - }); - expect(await finalHeightOfExperiment.textContent()).toMatch( - /My experiment 1: \d+ m/, - ); + // // Check that shared experiment has been run by + // // checking height in final height analysis + // const finalHeightAnalysis = sharedPage.getByRole("article", { + // name: "Final height", + // }); + // const finalHeightOfExperiment = finalHeightAnalysis.getByRole("listitem", { + // name: "My experiment 1", + // exact: true, + // }); + // expect(await finalHeightOfExperiment.textContent()).toMatch( + // /My experiment 1: \d+ m/, + // ); + // TODO: implement alternative check to see that experiment finished }); test("Given large app state, sharing is not possible", async ({ page }) => {