diff --git a/Sources/SpectroEDX/Fit/EnergyAxisRefinement.swift b/Sources/SpectroEDX/Fit/EnergyAxisRefinement.swift index ff3bc20..f07e61d 100644 --- a/Sources/SpectroEDX/Fit/EnergyAxisRefinement.swift +++ b/Sources/SpectroEDX/Fit/EnergyAxisRefinement.swift @@ -22,7 +22,17 @@ import Foundation import SpectroCore /// Which of the energy axis's offset, gain and the detector width the refinement adjusts, and how long it searches. Defaults are the -/// library's choice (all three, 600 simplex iterations, tolerance 1e-12, widths between 20 and 500 eV). +/// library's choice (all three, 600 simplex iterations, tolerance 1e-12, widths between 20 and 500 eV, offset within +-200 eV and +/// gain within +-5 % of the file's). +/// +/// THE OFFSET AND GAIN BOUNDS (0.13.2). A spectrum with one strong line fixes two numbers (its channel and its width in channels) +/// against three parameters, so the residual has a flat valley along which offset, gain and width trade; the simplex then walks +/// until a bound stops it (a synthetic one-line spectrum went to gain +444 %, offset -6.6 keV). +-200 eV and +-5 % are a stated +/// choice, far outside a working detector calibration (the tests plant up to +-12 eV and 0.2 %), and they keep the fitted +/// energies near or above zero from the 0.2 keV start. A refinement that ends against any bound is not a minimum: it is reported not +/// converged, with the parameter named, and the pooled fit uses the file's axis. +/// A policy stored before 0.13.2 has no offset or gain key and reads these defaults, not "no bound": unlike the rule for other new +/// `FitPolicy` keys, a re-run can differ from its stored run, and only where that run walked against the new bounds. public nonisolated struct AxisRefinementSettings: Sendable, Equatable, Codable { /// Adjust the offset. Default true; false keeps the file's. public var refineOffset = true @@ -38,12 +48,36 @@ public nonisolated struct AxisRefinementSettings: Sendable, Equatable, Codable { public var minimumResolutionEV = 20.0 /// A trial width at or above this, eV, is rejected. Default 500. public var maximumResolutionEV = 500.0 + /// A trial offset shift whose size is at or above this, eV, is rejected. Default 200. + public var maximumOffsetShiftEV = 200.0 + /// A trial relative gain shift whose size is at or above this is rejected. Default 0.05 (5 %). + public var maximumGainShift = 0.05 public init(refineOffset: Bool = true, refineGain: Bool = true, refineResolution: Bool = true, maxIterations: Int = 600, - tolerance: Double = 1e-12, minimumResolutionEV: Double = 20, maximumResolutionEV: Double = 500) { + tolerance: Double = 1e-12, minimumResolutionEV: Double = 20, maximumResolutionEV: Double = 500, + maximumOffsetShiftEV: Double = 200, maximumGainShift: Double = 0.05) { self.refineOffset = refineOffset; self.refineGain = refineGain; self.refineResolution = refineResolution self.maxIterations = maxIterations; self.tolerance = tolerance self.minimumResolutionEV = minimumResolutionEV; self.maximumResolutionEV = maximumResolutionEV + self.maximumOffsetShiftEV = maximumOffsetShiftEV; self.maximumGainShift = maximumGainShift + } + + private enum CodingKeys: String, CodingKey { + case refineOffset, refineGain, refineResolution, maxIterations, tolerance, minimumResolutionEV, maximumResolutionEV, + maximumOffsetShiftEV, maximumGainShift + } + + /// A policy stored before 0.13.2 has no offset or gain bound: it reads the defaults (see the type's note). + public init(from decoder: any Decoder) throws { + let c = try decoder.container(keyedBy: CodingKeys.self) + let defaults = AxisRefinementSettings() + self.init(refineOffset: try c.decode(Bool.self, forKey: .refineOffset), refineGain: try c.decode(Bool.self, forKey: .refineGain), + refineResolution: try c.decode(Bool.self, forKey: .refineResolution), + maxIterations: try c.decode(Int.self, forKey: .maxIterations), tolerance: try c.decode(Double.self, forKey: .tolerance), + minimumResolutionEV: try c.decode(Double.self, forKey: .minimumResolutionEV), + maximumResolutionEV: try c.decode(Double.self, forKey: .maximumResolutionEV), + maximumOffsetShiftEV: try c.decodeIfPresent(Double.self, forKey: .maximumOffsetShiftEV) ?? defaults.maximumOffsetShiftEV, + maximumGainShift: try c.decodeIfPresent(Double.self, forKey: .maximumGainShift) ?? defaults.maximumGainShift) } } @@ -58,7 +92,13 @@ public nonisolated struct AxisRefinementResult: Sendable, Equatable { package let rssFile: Double package let rssRefined: Double package let iterations: Int + /// The simplex met its tolerance AND no refined parameter ended against its bound. package let converged: Bool + /// The refined parameters that ended against their bound ("offset", "gain", "FWHM"), in that order; empty when none did. + /// Not empty: the result is not a minimum of the residual and the pooled fit does not use it. + public let parametersAtBound: [String] + /// "gain +5.000 % (bound +-5 %)", one clause per parameter at its bound; empty when none is. + package let boundText: String public var refinedAxis: EnergyAxis { EnergyAxis(offset: offset, scale: scale, size: size) } package let size: Int @@ -90,7 +130,8 @@ package nonisolated enum EnergyAxisRefinement { let d = (options.refineOffset ? u[0] : 0) * 1e-3 let g = (options.refineGain ? u[1] : 0) * 1e-4 let r = fileResolutionMnKaEV + (options.refineResolution ? u[2] : 0) - guard r > options.minimumResolutionEV, r < options.maximumResolutionEV else { return .infinity } + guard r > options.minimumResolutionEV, r < options.maximumResolutionEV, + abs(d * 1000) < options.maximumOffsetShiftEV, abs(g) < options.maximumGainShift else { return .infinity } let axis = EnergyAxis(offset: fileAxis.offset + d, scale: fileAxis.scale * (1 + g), size: fileAxis.size) let model = EDSLineModel.build(elements: settings.elements, axis: axis, beamEnergy: settings.beamEnergy, resolutionMnKaEV: r, escapePeaks: settings.escapePeaks, lineTail: settings.lineTail, @@ -117,12 +158,28 @@ package nonisolated enum EnergyAxisRefinement { let d = (options.refineOffset ? best[0] : 0) * 1e-3 let g = (options.refineGain ? best[1] : 0) * 1e-4 let r = fileResolutionMnKaEV + (options.refineResolution ? best[2] : 0) + // A parameter within 1e-4 of its bound span of a bound ended against it (the simplex presses on a rejected region). + var atBound: [String] = [], clauses: [String] = [] + func near(_ v: Double, _ lo: Double, _ hi: Double) -> Bool { min(v - lo, hi - v) <= 1e-4 * (hi - lo) } + if options.refineOffset, near(d * 1000, -options.maximumOffsetShiftEV, options.maximumOffsetShiftEV) { + atBound.append("offset") + clauses.append(String(format: "offset %+.1f eV (bound \u{00B1}%g eV)", d * 1000, options.maximumOffsetShiftEV)) + } + if options.refineGain, near(g, -options.maximumGainShift, options.maximumGainShift) { + atBound.append("gain") + clauses.append(String(format: "gain %+.3f %% (bound \u{00B1}%g %%)", g * 100, options.maximumGainShift * 100)) + } + if options.refineResolution, near(r, options.minimumResolutionEV, options.maximumResolutionEV) { + atBound.append("FWHM") + clauses.append(String(format: "FWHM %.0f eV (bounds %g-%g eV)", r, options.minimumResolutionEV, options.maximumResolutionEV)) + } try SpectroCancellation.check(cancelled) progress?(1) return AxisRefinementResult( fileOffset: fileAxis.offset, fileScale: fileAxis.scale, fileResolutionMnKaEV: fileResolutionMnKaEV, offset: fileAxis.offset + d, scale: fileAxis.scale * (1 + g), resolutionMnKaEV: r, - rssFile: rss([0, 0, 0]), rssRefined: rss(best), iterations: iterations, converged: converged, + rssFile: rss([0, 0, 0]), rssRefined: rss(best), iterations: iterations, converged: converged && atBound.isEmpty, + parametersAtBound: atBound, boundText: clauses.joined(separator: ", "), size: fileAxis.size) } diff --git a/Sources/SpectroEDX/Quant/PooledQuantification.swift b/Sources/SpectroEDX/Quant/PooledQuantification.swift index 519ba34..f659430 100644 --- a/Sources/SpectroEDX/Quant/PooledQuantification.swift +++ b/Sources/SpectroEDX/Quant/PooledQuantification.swift @@ -38,7 +38,9 @@ // WHAT IS NOT ESTABLISHED (read before trusting a number): // * The axis refinement was REFUTED at low counts: on demo-edx-tiny (706 counts, 24 px) it found -3.0 eV / +0.58 % // gain against the planted +10 eV / -0.2 %. The sigma terms carry NO axis term, so a wrong refined axis is not in -// the shown sigma. Its recovery against counts is unmeasured (open-items). +// the shown sigma. Its recovery against counts is unmeasured (open-items). A spectrum with one strong line does not +// determine offset, gain and width together (issue 83); the refinement then ends against a bound +// (`AxisRefinementSettings`), is reported not converged, and the file's axis is used, named in the warnings and footer. // * The absorption test (`SpectroscopyQuantifyTests`) checks plumbing only (applied/refused, terms present, the // composition moves); it does not check the correction against a truth geometry (open-items). // @@ -283,7 +285,13 @@ public nonisolated enum PooledQuantifier { catch { warnings.append("axis refinement failed (\(error)); the file's axis is used") } } if let r = refinement { - if r.rssRefined < r.rssFile { axis = r.refinedAxis; resolution = r.resolutionMnKaEV } + // A refinement that ended against a bound is not a minimum (issue 83: a one-line spectrum walked to gain +444 %). + if !r.parametersAtBound.isEmpty { + warnings.append("axis refinement stopped at its bound (\(r.boundText)): the spectrum does not determine the energy axis; the file's axis is used") + } else if r.rssRefined < r.rssFile { + axis = r.refinedAxis; resolution = r.resolutionMnKaEV + if !r.converged { warnings.append("axis refinement reached its iteration cap (\(r.iterations)) without converging; the refined axis is used") } + } else { warnings.append("axis refinement did not lower the residual; the file's axis is used") } } } @@ -509,7 +517,8 @@ public nonisolated enum PooledQuantifier { else if let r = refinement { footer.append(String(format: "energy axis: file %.3f keV + %.4f keV/ch \u{2192} refined %+.1f eV, gain %+.3f %%, FWHM(Mn K\u{03B1}) %.0f eV%@", r.fileOffset, r.fileScale, r.offsetShiftEV, r.gainShift * 100, r.resolutionMnKaEV, - axis == r.refinedAxis ? "" : "; NOT used (residual not lower)")) + !r.parametersAtBound.isEmpty ? "; NOT used (stopped at its bound: \(r.boundText))" + : axis == r.refinedAxis ? "" : "; NOT used (residual not lower)")) } if let b = BeamEnergy.provenance(typedKeV: input.method.beamEnergyKeV, metadata: input.metadata) { footer.append(b) } footer.append("\(fitRange.footerText) \u{00B7} \(input.pixelCount) px pooled (\(input.regionName))") diff --git a/Sources/SpectroEDX/Quant/ProposalAxisRefinement.swift b/Sources/SpectroEDX/Quant/ProposalAxisRefinement.swift index 1657b4b..34e4af2 100644 --- a/Sources/SpectroEDX/Quant/ProposalAxisRefinement.swift +++ b/Sources/SpectroEDX/Quant/ProposalAxisRefinement.swift @@ -80,7 +80,7 @@ extension ElementProposer { do { fit = try PooledQuantifier.run(input, tables: tables, cancelled: cancelled) } catch let c as SpectroCancelled { throw c } catch { break } - guard let r = fit.refinement, fit.usedAxis == r.refinedAxis else { break } + guard let r = fit.refinement, r.parametersAtBound.isEmpty, fit.usedAxis == r.refinedAxis else { break } anchorsUsed = anchors refinement = r axis = r.refinedAxis diff --git a/Sources/SpectroEDX/Quant/QuantificationPolicy.swift b/Sources/SpectroEDX/Quant/QuantificationPolicy.swift index 277dc47..58d617f 100644 --- a/Sources/SpectroEDX/Quant/QuantificationPolicy.swift +++ b/Sources/SpectroEDX/Quant/QuantificationPolicy.swift @@ -220,6 +220,10 @@ public nonisolated struct FitPolicy: Sendable, Equatable, Codable { || !(axisRefinement.minimumResolutionEV < axisRefinement.maximumResolutionEV) { return "The axis refinement needs a cap of 1 or more, a positive tolerance and a width range with its lower bound below its upper." } + if !(axisRefinement.maximumOffsetShiftEV.isFinite && axisRefinement.maximumOffsetShiftEV > 0) + || !(axisRefinement.maximumGainShift.isFinite && axisRefinement.maximumGainShift > 0) { + return "The axis refinement needs a finite, positive bound on the offset shift and on the gain shift." + } return nil } } diff --git a/Tests/SpectroCoreTests/Fixtures/spectro-golden/demo-edx-tiny-dm4.txt b/Tests/SpectroCoreTests/Fixtures/spectro-golden/demo-edx-tiny-dm4.txt index 4319c9e..14f7a6a 100644 --- a/Tests/SpectroCoreTests/Fixtures/spectro-golden/demo-edx-tiny-dm4.txt +++ b/Tests/SpectroCoreTests/Fixtures/spectro-golden/demo-edx-tiny-dm4.txt @@ -604,6 +604,8 @@ whole.quantification.refinement.rssFile = 40954ed514f96f63 1363.708087823318 whole.quantification.refinement.rssRefined = 40707f466a600e82 263.95469129106652 whole.quantification.refinement.iterations = 104 whole.quantification.refinement.converged = true +whole.quantification.refinement.parametersAtBound: n 0 sum 0000000000000000 0 fnv cbf29ce484222325 +whole.quantification.refinement.boundText = "" whole.quantification.refinement.size = 512 whole.quantification.plotModel: n 512 sum 408642460c1ad601 712.28420277562702 fnv 724e54731dbd18d7 whole.quantification.plotBackground: n 512 sum 406198f0ded13ffe 140.77940312260756 fnv dfb53193382aa4b diff --git a/Tests/SpectroCoreTests/Fixtures/spectro-golden/eds-unlisted-sim.txt b/Tests/SpectroCoreTests/Fixtures/spectro-golden/eds-unlisted-sim.txt index 19a675d..44b5b28 100644 --- a/Tests/SpectroCoreTests/Fixtures/spectro-golden/eds-unlisted-sim.txt +++ b/Tests/SpectroCoreTests/Fixtures/spectro-golden/eds-unlisted-sim.txt @@ -236,6 +236,8 @@ seed0.four.quantification.refinement.rssFile = 4224077700731f16 43012292665.5607 seed0.four.quantification.refinement.rssRefined = 42240444225b2cd6 42985460013.58757 seed0.four.quantification.refinement.iterations = 71 seed0.four.quantification.refinement.converged = true +seed0.four.quantification.refinement.parametersAtBound: n 0 sum 0000000000000000 0 fnv cbf29ce484222325 +seed0.four.quantification.refinement.boundText = "" seed0.four.quantification.refinement.size = 4096 seed0.four.quantification.plotModel: n 4096 sum 4150277a001e5418 4234728.0018511042 fnv 35bd38a650d42ab9 seed0.four.quantification.plotBackground: n 4096 sum 41443c04d862488b 2652169.6904993705 fnv 733b41c1d521b3d4 @@ -1483,6 +1485,8 @@ seed0.seven.quantification.refinement.rssFile = 4175186260e1f62a 22119974.055166 seed0.seven.quantification.refinement.rssRefined = 414a4bdcb7b8a233 3446713.43532207 seed0.seven.quantification.refinement.iterations = 94 seed0.seven.quantification.refinement.converged = true +seed0.seven.quantification.refinement.parametersAtBound: n 0 sum 0000000000000000 0 fnv cbf29ce484222325 +seed0.seven.quantification.refinement.boundText = "" seed0.seven.quantification.refinement.size = 4096 seed0.seven.quantification.plotModel: n 4096 sum 414e5cb35a013312 3979622.7031616056 fnv 19bc67e1afeb6e15 seed0.seven.quantification.plotBackground: n 4096 sum 412848dff77ef632 795759.9833905159 fnv fc0ed317795d391c diff --git a/Tests/SpectroCoreTests/Fixtures/spectro-golden/eds-weakline-sim.txt b/Tests/SpectroCoreTests/Fixtures/spectro-golden/eds-weakline-sim.txt index 9780841..7876d99 100644 --- a/Tests/SpectroCoreTests/Fixtures/spectro-golden/eds-weakline-sim.txt +++ b/Tests/SpectroCoreTests/Fixtures/spectro-golden/eds-weakline-sim.txt @@ -174,6 +174,8 @@ level0.quantification.refinement.rssFile = 414be7422632f3e3 3657348.2984299525 level0.quantification.refinement.rssRefined = 414b8fa4aa1546b0 3612489.3287742957 level0.quantification.refinement.iterations = 98 level0.quantification.refinement.converged = true +level0.quantification.refinement.parametersAtBound: n 0 sum 0000000000000000 0 fnv cbf29ce484222325 +level0.quantification.refinement.boundText = "" level0.quantification.refinement.size = 4096 level0.quantification.plotModel: n 4096 sum 4117b46554736a93 388377.33247152827 fnv 9185698d1dd3d97b level0.quantification.plotBackground: n 4096 sum 40f606b5aa5f6d78 90219.354094913113 fnv 6d86cfd7139ceb36 @@ -405,6 +407,8 @@ level4.quantification.refinement.rssFile = 413eb2e2eef45d2f 2011874.9334161987 level4.quantification.refinement.rssRefined = 413e42935bd6e312 1983123.3587476653 level4.quantification.refinement.iterations = 84 level4.quantification.refinement.converged = true +level4.quantification.refinement.parametersAtBound: n 0 sum 0000000000000000 0 fnv cbf29ce484222325 +level4.quantification.refinement.boundText = "" level4.quantification.refinement.size = 4096 level4.quantification.plotModel: n 4096 sum 4116312abc4c0005 363594.68388366728 fnv a01fd9b34b164cb9 level4.quantification.plotBackground: n 4096 sum 40f56682454b6097 87656.141917588524 fnv 5fd90c12681b255f diff --git a/Tests/SpectroEDXTests/AxisRefinementRunawayTests.swift b/Tests/SpectroEDXTests/AxisRefinementRunawayTests.swift new file mode 100644 index 0000000..f40d273 --- /dev/null +++ b/Tests/SpectroEDXTests/AxisRefinementRunawayTests.swift @@ -0,0 +1,121 @@ +// +// AxisRefinementRunawayTests.swift +// Role: The pooled energy-axis refinement on a synthetic spectrum with one strong line (Al K-alpha; Mg and Si K-alpha planted +// weak, both fitted at 0) on a 0-10.24 keV axis, 4322 counts in all. Under the broadened continuum edge step the +// refinement walked to gain +444 %, offset -6.6 keV and the 500 eV width bound and said converged; under the sharp step it +// stopped at gain +32.7 %, offset -484 eV. The file axis is the truth (the lines were planted on it). +// +// Cause: one line fixes two numbers (its channel and its width in channels), the refinement has three. Along the valley that +// keeps Al K-alpha at channel 148.6 and 9 channels wide the residual changes by under 5 % from gain 0 to +444 % (`valley`), so +// the step form only decides where on the valley the simplex stops. Refuting check run first: holding the continuum on the +// file axis and width left the broadened result bit-identical, so the continuum following the axis is not the cause. +// + +import Foundation +import Testing +import SpectroCore +@testable import SpectroEDX + +enum DemoSpectrum { + /// The pooled sum of a synthetic 32 x 32 x 1024 spectrum image (a deterministic LCG; offset 0, 10 eV/channel): a weak decaying + /// background and three Gaussian lines at 1.254, 1.487 and 1.740 keV (sigma^2 = 0.0035 keV^2), stronger inside a central + /// disc. Each pixel's channel value is truncated to an integer, then summed. Synthetic; the same generator as an app's demo + /// spectrum image. + static func pooled() -> [UInt64] { + let ny = 32, nx = 32, channels = 1024 + var sum = [UInt64](repeating: 0, count: channels) + var state: UInt64 = 0x9E3779B97F4A7C15 + func next() -> Double { + state = state &* 6364136223846793005 &+ 1442695040888963407 + return Double(state >> 40) / Double(1 << 24) + } + let peaks: [(keV: Double, matrix: Double, blob: Double)] = [(1.487, 40, 12), (1.254, 0.4, 14), (1.740, 0.4, 16)] + for y in 0.. PooledQuantificationInput { + var m = QuantificationMethod() + m.elements = ["Al", "Mg", "Si"].map { QuantificationMethod.ElementState(symbol: $0, role: .quantify, isManual: false) } + m.beamEnergyKeV = 200 + m.lockEnergyAxis = locked + let meta = SpectrumImageMetadata(fileName: "demo", filePath: "", scanWidth: 32, scanHeight: 32, channelCount: 1024, + energyOffsetEV: 0, energyDispersionEV: 10) + var i = PooledQuantificationInput(counts: pooled(), axis: axis, method: m, metadata: meta, regionName: "Whole map", + pixelCount: 1024) + i.policy = policy + return i + } +} + +struct AxisRefinementRunawayTests { + + private static func sharp() -> QuantificationPolicy { + var p = QuantificationPolicy.default; p.fit.continuumEdgeStep = .sharp; return p + } + + /// The fixture is the one the issue was reported on: 4322 counts. + @Test func fixtureIsTheReportedSpectrum() { + #expect(DemoSpectrum.pooled().reduce(0, +) == 4322) + } + + /// Under either edge step the spectrum does not determine the axis: the refinement stops at a bound, says so, is not + /// converged, and the fit runs on the file axis and width, number for number the locked-axis fit. + @Test(arguments: ["broadened", "sharp"]) + func oneLineSpectrumRefusesTheRefinedAxisWithANamedReason(_ step: String) throws { + let policy = step == "sharp" ? Self.sharp() : QuantificationPolicy.default + let q = try PooledQuantifier.run(DemoSpectrum.input(policy: policy), tables: nil) + let r = try #require(q.refinement) + #expect(!r.converged, "\(step): \(r.offsetShiftEV) eV, \(r.gainShift), \(r.resolutionMnKaEV) eV") + #expect(!r.parametersAtBound.isEmpty) + #expect(abs(r.gainShift) <= 0.05 + 1e-12 && abs(r.offsetShiftEV) <= 200 + 1e-9 && r.resolutionMnKaEV < 500) + #expect(q.usedAxis == DemoSpectrum.axis) + #expect(q.warnings.contains { $0.contains("stopped at its bound") && $0.contains("the file's axis is used") }, "\(q.warnings)") + #expect(q.footerLines.contains { $0.hasPrefix("energy axis:") && $0.contains("NOT used (stopped at its bound:") }, "\(q.footerLines)") + + let locked = try PooledQuantifier.run(DemoSpectrum.input(policy: policy, locked: true), tables: nil) + #expect(q.rows.map(\.net) == locked.rows.map(\.net)) + #expect(q.rows.map(\.sigma) == locked.rows.map(\.sigma)) + } + + /// The refinement objective at (d keV, g, R), on the file channels, as `EnergyAxisRefinement` computes it. + private static func rss(_ counts: [Double], step: ContinuumEdgeStep, d: Double, g: Double, r: Double) -> Double { + var fp = FitPolicy.default; fp.continuumEdgeStep = step + let file = DemoSpectrum.axis + let s = FitSettings.standard(elements: ["Al", "Mg", "Si"], axis: file, resolutionMnKaEV: 130, beamEnergy: 200, policy: fp) + let channels = file.fitChannels(from: s.fitFrom, to: s.fitTo) + let axis = EnergyAxis(offset: file.offset + d, scale: file.scale * (1 + g), size: file.size) + let model = EDSLineModel.build(elements: s.elements, axis: axis, beamEnergy: 200, resolutionMnKaEV: r, escapePeaks: s.escapePeaks) + let design = LinearDesign.build(model: model, axis: axis, channels: channels, background: s.background) + let y = channels.map { counts[$0] } + let sol = LeastSquaresFit.solve(nonNegative: design.nonNegative, free: design.free, y: y) + let mu = design.nonNegative.times(sol.nonNegative) + return y.indices.reduce(0) { $0 + (y[$1] - mu[$1]) * (y[$1] - mu[$1]) } + } + + /// The diagnosis: hold Al K-alpha (1.4865 keV) at channel 148.6 and its FWHM at 9.0 channels while the gain walks from 0 to + /// +444 %; the residual moves by under 5 % and, under the broadened step, falls all the way to the width bound. + @Test func oneLineValleyIsFlat() { + let counts = DemoSpectrum.pooled().map { Double($0) } + let gains = [0.0, 0.327, 1.5, 4.44] + var values: [Double] = [] + for g in gains { + let scale = 0.01 * (1 + g) + let fwhm = 9.0 * scale * 1000 + let r = (fwhm * fwhm - 2.5 * (1486.5 - 5898.7)).squareRoot() // the kit width law, FWHM^2 = R^2 + 2.5 (E - E_MnKa) eV^2 + values.append(Self.rss(counts, step: .broadened, d: 1.4865 - scale * 148.6, g: g, r: r)) + } + #expect((values.max()! - values.min()!) / values.max()! < 0.05, "\(values)") + #expect(zip(values, values.dropFirst()).allSatisfy { $0 > $1 }, "\(values)") + } +} diff --git a/Tests/SpectroEDXTests/QuantificationPolicyTests.swift b/Tests/SpectroEDXTests/QuantificationPolicyTests.swift index 87366fa..59b2c54 100644 --- a/Tests/SpectroEDXTests/QuantificationPolicyTests.swift +++ b/Tests/SpectroEDXTests/QuantificationPolicyTests.swift @@ -324,6 +324,21 @@ struct QuantificationPolicyTests { #expect(loose.iterations <= d.iterations) let bounded = try refined(AxisRefinementSettings(minimumResolutionEV: 20, maximumResolutionEV: 140)) #expect(bounded.resolutionMnKaEV < 140 && bounded.resolutionMnKaEV < d.resolutionMnKaEV, "\(bounded.resolutionMnKaEV) \(d.resolutionMnKaEV)") + // A refinement that ends against a bound is not converged, and names the parameter (issue 83). + #expect(d.converged && d.parametersAtBound.isEmpty) + #expect(!bounded.converged && bounded.parametersAtBound == ["FWHM"], "\(bounded.parametersAtBound)") + // The pooled fit does not use it: the file axis and width, and the footer says why. + let q = try Self.run(c, policy: Self.fit(FitPolicy(axisRefinement: AxisRefinementSettings(maximumResolutionEV: 140))), locked: false) + #expect(q.usedAxis == Self.axis) + #expect(q.footerLines.contains { $0.contains("NOT used (stopped at its bound: FWHM") }, "\(q.footerLines)") + // A capped refinement that is still moving is used, and says it did not converge. + let cappedRun = try Self.run(c, policy: Self.fit(FitPolicy(axisRefinement: AxisRefinementSettings(maxIterations: 3))), locked: false) + #expect(cappedRun.warnings.contains { $0.contains("iteration cap (3) without converging") }, "\(cappedRun.warnings)") + // The offset and gain bounds must be finite and positive. + for bad in [AxisRefinementSettings(maximumOffsetShiftEV: 0), AxisRefinementSettings(maximumGainShift: -0.01), + AxisRefinementSettings(maximumGainShift: .infinity)] { + #expect(FitPolicy(axisRefinement: bad).problem != nil) + } } // MARK: F13 the range-sensitivity alternative