diff --git a/src/__tests__/__snapshots__/index.test.ts.snap b/src/__tests__/__snapshots__/index.test.ts.snap index 71f924ef..055e9cf5 100644 --- a/src/__tests__/__snapshots__/index.test.ts.snap +++ b/src/__tests__/__snapshots__/index.test.ts.snap @@ -113,6 +113,7 @@ exports[`test existence of exported functions 1`] = ` "xyRealMinYPoint", "xyReduceNonContinuous", "xyRolling", + "xyRollingCircleTransform", "xySetYValue", "xySortX", "xyToXYArray", diff --git a/src/xy/__tests__/xyRollingCircleTransform.test.ts b/src/xy/__tests__/xyRollingCircleTransform.test.ts new file mode 100644 index 00000000..57559ed0 --- /dev/null +++ b/src/xy/__tests__/xyRollingCircleTransform.test.ts @@ -0,0 +1,99 @@ +import type { DataXY } from 'cheminfo-types'; +import { expect, test } from 'vitest'; + +import { xyRollingCircleTransform } from '../xyRollingCircleTransform.ts'; + +test('simple slope', () => { + const data: DataXY = { + x: [], + y: [], + }; + for (let i = 0; i < 5; i++) { + data.x.push(i); + data.y.push(i); + } + + const result1 = xyRollingCircleTransform(data); + + expect(result1).toStrictEqual(Float64Array.from([1, 2, 3, 4, 5])); + + const result2 = xyRollingCircleTransform(data, { radius: 2 }); + + expect(result2).toStrictEqual( + Float64Array.from([ + 2.732050807568877, 3.732050807568877, 4.732050807568877, + 5.732050807568877, 6, + ]), + ); + + const result3 = xyRollingCircleTransform(data, { + radius: 1, + shifted: false, + }); + + expect(result3).toStrictEqual(Float64Array.from([0, 1, 2, 3, 4])); + + const result4 = xyRollingCircleTransform(data, { + radius: 1, + position: 'bottom', + }); + + expect(result4).toStrictEqual(Float64Array.from([-1, -0, 1, 2, 3])); + + const result5 = xyRollingCircleTransform(data, { + position: 'bottom', + shifted: false, + }); + + expect(result5).toStrictEqual(Float64Array.from([-0, 1, 2, 3, 4])); +}); + +test('steep slope', () => { + const data: DataXY = { + x: [], + y: [], + }; + for (let i = 0; i < 5; i++) { + data.x.push(i); + data.y.push(i * 10); + } + const result1 = xyRollingCircleTransform(data, { + shifted: false, + radius: 2, + }); + + expect(result1).toStrictEqual( + Float64Array.from([18, 28, 38, 39.732050807568875, 40]), + ); +}); + +test('wrong position', () => { + const data: DataXY = { + x: [0, 1, 2], + y: [0, 1, 2], + }; + + expect(() => { + xyRollingCircleTransform(data, { position: 'middle' as 'top' }); + }).toThrow('Invalid position: middle'); +}); + +test('two peaks', () => { + const data: DataXY = { + x: [0, 1, 2, 3, 4, 5, 6, 7, 8], + y: [0, 1, 2, 1, 0, 1, 2, 1, 0], + }; + + const result = xyRollingCircleTransform(data, { + shifted: false, + radius: 5, + }); + + expect(result).toStrictEqual( + Float64Array.from([ + 1.5825756949558398, 1.8989794855663558, 2, 1.8989794855663558, + 1.5825756949558398, 1.8989794855663558, 2, 1.8989794855663558, + 1.5825756949558398, + ]), + ); +}); diff --git a/src/xy/index.ts b/src/xy/index.ts index 6f8ecc1f..78d4bcff 100644 --- a/src/xy/index.ts +++ b/src/xy/index.ts @@ -35,6 +35,7 @@ export * from './xyRealMinYPoint.ts'; export { xyReduce } from './xyReduce.ts'; export * from './xyReduceNonContinuous.ts'; export * from './xyRolling.ts'; +export * from './xyRollingCircleTransform.ts'; export * from './xySetYValue.ts'; export * from './xySortX.ts'; export * from './xyToXYArray.ts'; diff --git a/src/xy/xyRollingCircleTransform.ts b/src/xy/xyRollingCircleTransform.ts new file mode 100644 index 00000000..5f384a33 --- /dev/null +++ b/src/xy/xyRollingCircleTransform.ts @@ -0,0 +1,91 @@ +/** + * + * Dong, Jian, et al. "An algorithm of filtering noises in multi-beam data based on rolling circle transform." 2019 2nd International Conference on Sustainable Energy, Environment and Information Engineering (SEEIE 2019). Atlantis Press, 2019. + * DONG Jian, PENG Rencan, ZHANG Lihua, WANG Zhijun. An Algorithm of Filtering Noises in Multi-beam Data Based on Rolling Circle Transform[J]. Geomatics and Information Science of Wuhan University, 2016, 41(1): 86-92. DOI: 10.13203/j.whugis20130757 + * @param data + * @param options + */ + +import type { DataXY } from 'cheminfo-types'; + +import { xFindClosestIndex } from '../x/xFindClosestIndex.ts'; + +interface XYRollingCircleTransformOptions { + /** + * The radius of the circle used for the rolling circle transform. + * It should be a positive number and in the 'x' unit + * @default 1 + */ + radius?: number; // radius of the circle + /** + * Should the circle scan the top or the bottom of the XY data ? + * @default 'top' + */ + position?: 'top' | 'bottom'; // position of the circle relative to the peak + /** + * Should we keep the Y centers of the circles or should we move the centers + * so that it touches the XY data ? + * @default: true + */ + shifted?: boolean; // if true, the y values will be shifted to fit the circle center +} + +export function xyRollingCircleTransform( + data: DataXY, + options: XYRollingCircleTransformOptions = {}, +): Float64Array { + const { x } = data; + let { y } = data; + const { radius = 1, position = 'top', shifted = true } = options; + + if (position !== 'top' && position !== 'bottom') { + throw new Error(`Invalid position: ${String(position)}`); + } + + if (position === 'bottom') { + y = y.slice(); + for (let i = 0; i < y.length; i++) { + y[i] = -y[i]; + } + } + + if (x.length === 0 || y.length === 0) { + return new Float64Array(); + } + const yCenters = new Float64Array(x.length); + for (let i = 0; i < x.length; i++) { + const x0 = x[i]; // x center of the current circle + const fromX = xFindClosestIndex(x, x0 - radius); + const toX = xFindClosestIndex(x, x0 + radius); + + // for the circle radius we need to evaluate the minimal vertical shift + + const y0 = y[i] + radius; // y center of the circle on top of peak + let yShift = y0; // this is the minimal possible shift + for (let j = fromX; j <= toX; j++) { + const currentX = x[j]; + if (currentX < x0 - radius || currentX > x0 + radius) { + continue; + } + const currentMinYShift = + y[j] + Math.sqrt(radius ** 2 - (currentX - x0) ** 2); + if (currentMinYShift > yShift) { + yShift = currentMinYShift; + } + } + yCenters[i] = yShift; + } + if (!shifted) { + for (let i = 0; i < yCenters.length; i++) { + yCenters[i] -= radius; + } + } + + if (position === 'bottom') { + for (let i = 0; i < yCenters.length; i++) { + yCenters[i] = -yCenters[i]; + } + } + + return yCenters; +}