{"spec_id":"stereonet-equal-area","library":"echarts","language":"javascript","code":"// anyplot.ai\n// stereonet-equal-area: Structural Geology Stereonet (Equal-Area Projection)\n// Library: echarts 5.5.1 | JavaScript 22.22.3\n// Quality: 90/100 | Created: 2026-06-16\n//# anyplot-orientation: square\n// anyplot.ai\n// stereonet-equal-area: Structural Geology Stereonet (Equal-Area Projection)\n// Library: echarts 5.5.1 | JavaScript\n// Quality: pending | Created: 2026-06-16\n//\n// Schmidt lower-hemisphere equal-area net rendered idiomatically in ECharts.\n// ECharts has no native stereonet, so we build the projection ourselves on a\n// hidden square cartesian2d grid: the equal-area net (meridians + small circles)\n// and the great circles / poles are computed in data coords [-1,1] and drawn with\n// `custom` series + `scatter`. Imprint palette throughout; only chrome flips with\n// the theme.\n\nconst t = window.ANYPLOT_TOKENS;\nconst THEME = window.ANYPLOT_THEME;\nconst ink = t.ink;\nconst inkSoft = t.inkSoft;\nconst grid = t.grid;\nconst muted = THEME === \"light\" ? \"#6B6A63\" : \"#A8A79F\";\n\nconst DEG = Math.PI / 180;\nconst SQRT2 = Math.SQRT2;\n\n// --- Deterministic RNG (LCG + Box-Muller) ----------------------------------\nlet _seed = 20260616 >>> 0;\nfunction rand() {\n  _seed = (_seed * 1664525 + 1013904223) >>> 0;\n  return _seed / 4294967296;\n}\nfunction gauss() {\n  let u = 0, v = 0;\n  while (u === 0) u = rand();\n  while (v === 0) v = rand();\n  return Math.sqrt(-2 * Math.log(u)) * Math.cos(2 * Math.PI * v);\n}\nconst clamp = (x, lo, hi) => Math.min(hi, Math.max(lo, x));\n\n// --- Equal-area (Schmidt) projection helpers -------------------------------\n// Convention: vectors are [E(east), N(north), D(down)]; lower hemisphere is D>=0.\n// The net centre is the vertical line (plunge 90); the primitive circle (r=1) is\n// the horizontal plane. Azimuth (trend) is measured clockwise from North (top).\nfunction vecToXY(v) {\n  let [E, N, D] = v;\n  if (D < 0) { E = -E; N = -N; D = -D; } // force lower hemisphere\n  const horiz = Math.hypot(E, N);\n  const theta = Math.acos(clamp(D, -1, 1)); // angular distance from centre\n  const r = SQRT2 * Math.sin(theta / 2);    // equal-area radius (=1 at horizon)\n  if (horiz < 1e-9) return [0, 0];\n  return [r * (E / horiz), r * (N / horiz)];\n}\n\n// Pole (normal) to a plane given strike & dip (right-hand rule: dip dir = strike+90)\nfunction poleVec(strike, dip) {\n  const dipDir = strike + 90;\n  const pTrend = (dipDir + 180) * DEG; // pole azimuth\n  const pPlunge = (90 - dip) * DEG;    // pole plunge below horizontal\n  const cp = Math.cos(pPlunge);\n  return [cp * Math.sin(pTrend), cp * Math.cos(pTrend), Math.sin(pPlunge)];\n}\n\nfunction norm(v) {\n  const m = Math.hypot(v[0], v[1], v[2]) || 1;\n  return [v[0] / m, v[1] / m, v[2] / m];\n}\nfunction cross(a, b) {\n  return [a[1] * b[2] - a[2] * b[1], a[2] * b[0] - a[0] * b[2], a[0] * b[1] - a[1] * b[0]];\n}\n\n// Great-circle arc of the plane whose pole is `n` (lower-hemisphere half).\nfunction greatCircleArc(n) {\n  n = norm(n);\n  const u = norm(cross(n, [0, 0, 1])); // horizontal in-plane direction\n  const v = cross(n, u);               // completes the in-plane basis\n  const pts = [];\n  for (let k = 0; k <= 180; k++) {\n    const phi = (-Math.PI / 2) + (k / 180) * Math.PI; // half-circle, contiguous\n    const D = Math.sin(phi) * v[2];\n    if (D < -1e-9) continue;\n    pts.push(vecToXY([\n      Math.cos(phi) * u[0] + Math.sin(phi) * v[0],\n      Math.cos(phi) * u[1] + Math.sin(phi) * v[1],\n      Math.cos(phi) * u[2] + Math.sin(phi) * v[2],\n    ]));\n  }\n  return pts;\n}\n\n// --- The equal-area net (subtle background grid) ---------------------------\nconst netLines = [];\n// Meridians: great circles sharing the N-S axis, dipping E and W every 10 deg.\nfor (let d = 10; d <= 80; d += 10) {\n  netLines.push(greatCircleArc(poleVec(0, d)));   // east-dipping\n  netLines.push(greatCircleArc(poleVec(180, d))); // west-dipping\n}\nnetLines.push(greatCircleArc(poleVec(0, 90)));    // N-S diameter\n// Small circles: cones about the horizontal N-S axis every 10 deg.\nconst nAxis = [0, 1, 0], e1 = [1, 0, 0], e2 = [0, 0, 1];\nfor (let a = 10; a <= 170; a += 10) {\n  const ca = Math.cos(a * DEG), sa = Math.sin(a * DEG);\n  const arc = [];\n  for (let k = 0; k <= 180; k++) {\n    const tt = (k / 180) * Math.PI;\n    const vec = [\n      ca * nAxis[0] + sa * (Math.cos(tt) * e1[0] + Math.sin(tt) * e2[0]),\n      ca * nAxis[1] + sa * (Math.cos(tt) * e1[1] + Math.sin(tt) * e2[1]),\n      ca * nAxis[2] + sa * (Math.cos(tt) * e1[2] + Math.sin(tt) * e2[2]),\n    ];\n    if (vec[2] < -1e-9) continue;\n    arc.push(vecToXY(vec));\n  }\n  if (arc.length > 1) netLines.push(arc);\n}\n// Primitive circle (horizon) + perimeter ticks + cardinal labels.\nconst primitive = [];\nfor (let k = 0; k <= 360; k++) primitive.push([Math.sin(k * DEG), Math.cos(k * DEG)]);\nconst ticks = [];\nfor (let az = 0; az < 360; az += 10) {\n  const outer = az % 90 === 0 ? 1.05 : az % 30 === 0 ? 1.035 : 1.022;\n  const s = Math.sin(az * DEG), c = Math.cos(az * DEG);\n  ticks.push([[s, c], [outer * s, outer * c]]);\n}\nconst cardinals = [\n  { t: \"N\", a: 0 }, { t: \"E\", a: 90 }, { t: \"S\", a: 180 }, { t: \"W\", a: 270 },\n];\n\n// --- Field data: poles to planes, by feature type --------------------------\n// Imprint palette, canonical order; faults take the matte-red semantic anchor.\nconst SETS = [\n  { name: \"Bedding\", strike: 42, dip: 24, sS: 11, sD: 7, n: 38, color: t.palette[0] },\n  { name: \"Joint set 1\", strike: 118, dip: 78, sS: 9, sD: 8, n: 32, color: t.palette[1] },\n  { name: \"Joint set 2\", strike: 205, dip: 66, sS: 12, sD: 9, n: 27, color: t.palette[2] },\n  { name: \"Fault\", strike: 302, dip: 52, sS: 14, sD: 10, n: 14, color: t.palette[4] },\n];\n\nconst allPoles = []; // projected [x,y] for density\nconst poleSeries = SETS.map((set) => {\n  const data = [];\n  let mean = [0, 0, 0];\n  for (let i = 0; i < set.n; i++) {\n    const strike = set.strike + gauss() * set.sS;\n    const dip = clamp(set.dip + gauss() * set.sD, 2, 88);\n    const pv = poleVec(strike, dip);\n    mean = [mean[0] + pv[0], mean[1] + pv[1], mean[2] + (pv[2] < 0 ? -pv[2] : pv[2])];\n    const xy = vecToXY(pv);\n    data.push(xy);\n    allPoles.push(xy);\n  }\n  return { ...set, data, mean: norm(mean) };\n});\n\n// --- Kamb-style density (Gaussian KDE in the projection plane) -------------\nconst N = 110, H = 0.13, INV = 1 / (2 * H * H);\nconst gx = [], gy = [];\nfor (let i = 0; i <= N; i++) {\n  const c = -1 + (2 * i) / N;\n  gx.push(c); gy.push(c);\n}\nconst Z = [];\nlet maxD = 0;\nfor (let j = 0; j <= N; j++) {\n  Z[j] = [];\n  for (let i = 0; i <= N; i++) {\n    const x = gx[i], y = gy[j];\n    if (x * x + y * y > 1.0) { Z[j][i] = NaN; continue; }\n    let s = 0;\n    for (const p of allPoles) {\n      const dx = x - p[0], dy = y - p[1];\n      s += Math.exp(-(dx * dx + dy * dy) * INV);\n    }\n    Z[j][i] = s;\n    if (s > maxD) maxD = s;\n  }\n}\n// Marching squares -> contour segments coloured along imprint_seq (low->high).\nfunction hex(c) { const n = parseInt(c.slice(1), 16); return [n >> 16 & 255, n >> 8 & 255, n & 255]; }\nconst seqLo = hex(t.seq[0]), seqHi = hex(t.seq[1]);\nfunction seqColor(f) {\n  const r = Math.round(seqLo[0] + (seqHi[0] - seqLo[0]) * f);\n  const g = Math.round(seqLo[1] + (seqHi[1] - seqLo[1]) * f);\n  const b = Math.round(seqLo[2] + (seqHi[2] - seqLo[2]) * f);\n  return `rgb(${r},${g},${b})`;\n}\nconst LEVF = [0.12, 0.25, 0.4, 0.55, 0.7, 0.85];\nconst contours = [];\nfor (let li = 0; li < LEVF.length; li++) {\n  const L = LEVF[li] * maxD;\n  const col = seqColor(li / (LEVF.length - 1));\n  const lw = 1.4 + 1.4 * (li / (LEVF.length - 1));\n  for (let j = 0; j < N; j++) {\n    for (let i = 0; i < N; i++) {\n      const v0 = Z[j][i], v1 = Z[j][i + 1], v2 = Z[j + 1][i + 1], v3 = Z[j + 1][i];\n      if (isNaN(v0) || isNaN(v1) || isNaN(v2) || isNaN(v3)) continue;\n      const pts = [];\n      if ((v0 >= L) !== (v1 >= L)) { const tA = (L - v0) / (v1 - v0); pts.push([gx[i] + tA * (gx[i + 1] - gx[i]), gy[j]]); }\n      if ((v1 >= L) !== (v2 >= L)) { const tA = (L - v1) / (v2 - v1); pts.push([gx[i + 1], gy[j] + tA * (gy[j + 1] - gy[j])]); }\n      if ((v3 >= L) !== (v2 >= L)) { const tA = (L - v3) / (v2 - v3); pts.push([gx[i] + tA * (gx[i + 1] - gx[i]), gy[j + 1]]); }\n      if ((v0 >= L) !== (v3 >= L)) { const tA = (L - v0) / (v3 - v0); pts.push([gx[i], gy[j] + tA * (gy[j + 1] - gy[j])]); }\n      if (pts.length === 2) contours.push({ p1: pts[0], p2: pts[1], col, lw });\n      else if (pts.length === 4) {\n        contours.push({ p1: pts[0], p2: pts[1], col, lw });\n        contours.push({ p1: pts[2], p2: pts[3], col, lw });\n      }\n    }\n  }\n}\n\n// --- Mean great circle per feature type (planar orientation) ---------------\nconst meanCircles = poleSeries.map((s) => ({ points: greatCircleArc(s.mean), color: s.color }));\n\n// --- Init & render ----------------------------------------------------------\nconst chart = echarts.init(document.getElementById(\"container\"));\nchart.on(\"finished\", () => { window.__anyplotReady = true; });\n\nconst AX = { type: \"value\", min: -1.18, max: 1.18, show: false };\n\n// custom-series renderers (silent background layers)\nfunction polylineChildren(lines, styleFn) {\n  return (params, api) => ({\n    type: \"group\",\n    silent: true,\n    children: lines.map((ln, idx) => ({\n      type: \"polyline\",\n      shape: { points: ln.points.map((p) => api.coord(p)) },\n      style: { stroke: ln.color, lineWidth: ln.lw, fill: \"none\", lineCap: \"round\" },\n      ...(styleFn ? styleFn(idx) : {}),\n    })),\n  });\n}\n\nconst netRender = (params, api) => {\n  const children = [];\n  // primitive circle\n  children.push({ type: \"polyline\", silent: true,\n    shape: { points: primitive.map((p) => api.coord(p)) },\n    style: { stroke: inkSoft, lineWidth: 2.4, fill: \"none\" } });\n  // net meridians + small circles\n  for (const ln of netLines) {\n    children.push({ type: \"polyline\", silent: true,\n      shape: { points: ln.map((p) => api.coord(p)) },\n      style: { stroke: grid, lineWidth: 1, fill: \"none\" } });\n  }\n  // perimeter ticks\n  for (const tk of ticks) {\n    children.push({ type: \"line\", silent: true,\n      shape: { x1: api.coord(tk[0])[0], y1: api.coord(tk[0])[1], x2: api.coord(tk[1])[0], y2: api.coord(tk[1])[1] },\n      style: { stroke: inkSoft, lineWidth: 1.6 } });\n  }\n  // cardinal labels\n  for (const cd of cardinals) {\n    const pos = api.coord([1.12 * Math.sin(cd.a * DEG), 1.12 * Math.cos(cd.a * DEG)]);\n    children.push({ type: \"text\", silent: true,\n      style: { text: cd.t, x: pos[0], y: pos[1], fill: cd.t === \"N\" ? ink : inkSoft,\n        fontSize: cd.t === \"N\" ? 26 : 20, fontWeight: cd.t === \"N\" ? \"bold\" : \"normal\",\n        align: \"center\", verticalAlign: \"middle\", fontFamily: \"sans-serif\" } });\n  }\n  // North arrow (small triangle above N)\n  const nb = api.coord([0, 1.04]), nt = api.coord([0, 1.10]);\n  children.push({ type: \"polygon\", silent: true,\n    shape: { points: [[nt[0], nt[1]], [nb[0] - 8, nb[1]], [nb[0] + 8, nb[1]]] },\n    style: { fill: ink } });\n  return { type: \"group\", children };\n};\n\nconst option = {\n  animation: false,\n  backgroundColor: \"transparent\",\n  color: t.palette,\n  title: {\n    text: \"stereonet-equal-area · javascript · echarts · anyplot.ai\",\n    subtext: \"Lower-hemisphere equal-area (Schmidt) net · poles to planes, Kamb density & mean great circles\",\n    left: \"center\", top: 18,\n    textStyle: { color: ink, fontSize: 22, fontWeight: \"bold\" },\n    subtextStyle: { color: inkSoft, fontSize: 15 },\n  },\n  legend: {\n    data: SETS.map((s) => s.name),\n    orient: \"vertical\", left: 36, top: 120,\n    itemWidth: 20, itemHeight: 14, itemGap: 14,\n    icon: \"circle\",\n    textStyle: { color: ink, fontSize: 17 },\n  },\n  grid: { left: 140, right: 140, top: 160, bottom: 120 },\n  xAxis: AX,\n  yAxis: AX,\n  series: [\n    // 1. equal-area net (background)\n    { type: \"custom\", coordinateSystem: \"cartesian2d\", silent: true, z: 1,\n      data: [0], renderItem: netRender },\n    // 2. Kamb density contours (imprint_seq, low->high)\n    { type: \"custom\", coordinateSystem: \"cartesian2d\", silent: true, z: 2,\n      data: [0],\n      renderItem: (params, api) => ({\n        type: \"group\", silent: true,\n        children: contours.map((c) => {\n          const a = api.coord(c.p1), b = api.coord(c.p2);\n          return { type: \"line\", silent: true,\n            shape: { x1: a[0], y1: a[1], x2: b[0], y2: b[1] },\n            style: { stroke: c.col, lineWidth: c.lw, opacity: 0.85 } };\n        }),\n      }) },\n    // 3. mean great circles per feature type\n    { type: \"custom\", coordinateSystem: \"cartesian2d\", silent: true, z: 3,\n      data: [0],\n      renderItem: (params, api) => ({\n        type: \"group\", silent: true,\n        children: meanCircles.map((mc) => ({\n          type: \"polyline\",\n          shape: { points: mc.points.map((p) => api.coord(p)) },\n          style: { stroke: mc.color, lineWidth: 3, fill: \"none\", opacity: 0.9, lineCap: \"round\" },\n        })),\n      }) },\n    // 4. poles to planes (points), one series per feature type -> legend\n    ...poleSeries.map((s) => ({\n      name: s.name, type: \"scatter\", coordinateSystem: \"cartesian2d\", z: 4,\n      data: s.data, symbolSize: 14,\n      itemStyle: { color: s.color, borderColor: t.pageBg, borderWidth: 1.2, opacity: 0.95 },\n    })),\n  ],\n};\n\nchart.setOption(option);\n"}