Engineering maths

Overview

The spectrum analyzer, the control chart and the strip chart display results that are also useful without the display: a vibration monitor computes a spectrum in a worker, a quality system computes Cpk on the server, and a historian client reduces a day of samples before it plots them. For this reason the calculations are provided in three modules of plain functions with no DOM dependency, and the components call these functions.

All three modules are tested against reference values: the Bluestein FFT against a direct DFT at non-power-of-two lengths, the window gains against their published values, and the SPC constants against the published tables. Each module is included in the bundle and is also available as its own file under source/modules/.

In a Worker and in Node

In a page, the files register the namespaces above on Smart; smart.dsp.js and smart.decimate.js do so when the bundle or another module has brought the framework first, and are otherwise available as window.SmartDSP and window.SmartDecimate. Where there is no document, in a Worker or in Node, the files run the functions alone, without the components and the framework, on a Smart object of their own:

//A classic Worker: importScripts() takes the files themselves.
importScripts('smart-industrial/source/modules/smart.dsp.js');        //self.Smart.DSP
importScripts('smart-industrial/source/modules/smart.decimate.js');   //self.Smart.Utilities.Decimate
importScripts('smart-industrial/source/modules/smart.spc.js');        //self.Smart.Industrial.spc

//Node: require() returns that Smart object.
const { Industrial } = require('smart-industrial/source/modules/smart.spc.js');
Industrial.spc.capability(values, { lsl: 495, usl: 505 });

//The bundle carries all three, and the Units and Audit modules.
const all = require('smart-industrial');   //all.DSP, all.Industrial.spc, all.Industrial.Units, ...

The paths are relative to wherever the application serves the package; a Worker resolves them against the Worker script's own URL.

Signal processing

const DSP = Smart.DSP;

const result = DSP.spectrum(samples, { sampleRate: 48000, window: 'hann', scaling: 'amplitude' });
//{ frequencies, magnitudes, binWidth, size, window } - one-sided, window gain corrected

const peaks = DSP.peaks(result, { count: 5, minProminence: 0 });
//[{ frequency, magnitude, index }], strongest first

const distortion = DSP.thd(samples, { sampleRate: 48000, fundamental: 1000, harmonics: 5 });
//{ thd, thdDb, fundamental, amplitudes } - thd is a ratio

const noiseFloor = DSP.psd(samples, { sampleRate: 48000, size: 1024, overlap: 0.5 });   //Welch
const frames = DSP.spectrogram(samples, { sampleRate: 48000, size: 512, overlap: 0.5 });
//{ frames, frequencies, times, binWidth }

DSP.rms(samples);
DSP.toDb(result.magnitudes, { reference: 1, floor: -200, kind: 'amplitude' });
DSP.window('blackman-harris', 4096);          //the coefficients
DSP.windowGain('hann');                        //{ coherent, noise } - both, on purpose
DSP.windows;                                   //['rectangular', 'hann', 'hamming', 'blackman', 'blackman-harris', 'flat-top']

//The transforms themselves, on separate real and imaginary arrays.
DSP.fft(re, im);            //in place, power-of-two length only (other lengths throw)
DSP.bluestein(re, im);      //in place, any length: the FFT when it can, Bluestein when it must
DSP.ifft(re, im);           //in place, scaled by 1/n

//A real signal, windowed and transformed. size defaults to the next power of two, so a
//signal of another length is zero-padded; pass size to transform exactly that many samples.
DSP.transform(samples, { window: 'hann', size: samples.length });   //{ re, im, window, size }

The scaling is explicit because the two window gains are different. A Hann-windowed sine reads half its true amplitude unless the result is divided by the coherent gain of the window, and a power density needs the noise gain instead. spectrum() takes a scaling argument of amplitude, power or density and applies the correct gain.

The functions do not use the DOM, so for transform sizes that would block the main thread the module can be loaded into a Worker, as described above.

Statistical process control

const spc = Smart.Industrial.spc;

//Variables charts. Each returns { type, points, indexes, center, ucl, lcl, sigma, sigmaWithin, n, secondary },
//secondary being the matching R, S or MR chart. indexes[i] is the sample points[i] came from: a sample
//with no usable value is left out of points, so the two can differ.
const chart = spc.xbarR(subgroups, { baselineCount: 20 });
spc.xbarS(subgroups);
spc.iMR(values);

//Attribute charts. With unequal sample sizes ucl and lcl are arrays, one per point.
spc.pChart(defectives, sizes);
spc.npChart(defectives, size);
spc.cChart(defects);
spc.uChart(defects, sizes);

//Which points break which rules: [{ index, point, rules: ['beyond-3-sigma', 'two-of-three-beyond-2-sigma'] }]
//point is the position in chart.points; index is the sample it came from, read from the indexes the
//chart is passed with. Pass the chart (or { center, sigma, indexes: chart.indexes }): without indexes,
//index is the same as point, and after a skipped sample it names the wrong sample.
const violations = spc.rules(chart.points, chart, 'nelson');
spc.ruleNames.nelson;           //the eight Nelson rule names; .westernElectric the four

//Capability. Cp and Cpk from sigma within; Pp and Ppk from sigma overall.
const study = spc.capability(subgroups, { lsl: 495, usl: 505 });
//{ n, mean, sigmaWithin, sigmaOverall, cp, cpk, pp, ppk, ppm }
//Subgroups (an array of arrays) take sigma within from the R chart, or from the S chart with
//{ type: 'xbar-s' }; a flat list of values takes it from the moving range.

spc.constants.d2(5);            //d2, d3, a2 and c4 as functions of the subgroup size
spc.stats.mean(values);  spc.stats.stdev(values);

• Cp and Pp are different measurements. Cp and Cpk use the within-subgroup sigma (the short-term spread, estimated from the subgroup ranges or, for single values, from the moving range); Pp and Ppk use the overall sigma. With only one specification limit there is no Cp, because a one-sided tolerance has no width, and the function returns null instead of zero.
• The limits of a p-chart depend on the sample size. A smaller sample gives less information and therefore wider limits. The limits are computed for each sample and drawn as steps.
• Rules are identified by name. The run rules are named, for example fourteen-alternating, instead of numbered, so a rule name is readable without the reference table.

Series reduction

const Decimate = Smart.Utilities.Decimate;

Decimate.minMax(points, plotWidth, { x: 'time', y: 'value' });     //[{ x, y }], at most two per bucket
Decimate.minMaxIndexed(count, function (index) { return ring[index]; }, plotWidth);   //indices, in place
Decimate.lttb(points, 500);                                        //[{ x, y }], 500 kept incl. endpoints
Decimate.stepPreserving(transitions, pixels, { from, to });        //[{ time, value, collapsed }]

Min/max is the default for waveforms because the envelope is the signal at one sample per pixel: a one-sample spike in fifty thousand samples still reaches the top of the chart. LTTB draws a smoother line and is suitable for a trend of averages, but it keeps one sample per bucket, so a peak and a trough a few samples apart are reduced to one of them. The step-preserving reduction is for digital signals, where a dropped edge would hide a state. The strip chart uses minMaxIndexed over its ring buffer.

Missing values

The chart and capability functions skip null, undefined, '', booleans and NaN instead of letting Number(null) convert them to 0. Otherwise a gap in an acquisition would become a measured zero: it would lower the mean, draw a column at the bottom of the scale and report aria-valuenow="0". The guard is in the modules so that every caller gets it. Two helpers are the exception: spc.stats.mean() and spc.stats.stdev() add what they are given and expect numbers only. The reduction functions treat the same values as gaps: minMaxIndexed() keeps the index of every gap in its result, which is how the strip chart lifts the pen there, while minMax() and lttb() leave gaps out, so a line drawn through their result joins the samples on either side. The DSP functions treat samples the same way: rms() leaves a gap out of the sum, and the windowed transform treats it as silence.