dispersive lorentzian implementation

src/index.ts

@@ -1,5 +1,6 @@
 export * from './shapes/1d/gaussian/Gaussian';
 export * from './shapes/1d/lorentzian/Lorentzian';
+export * from './shapes/1d/lorentzianDispersive/LorentzianDispersive';
 export * from './shapes/1d/pseudoVoigt/PseudoVoigt';
 export * from './shapes/1d/generalizedLorentzian/GeneralizedLorentzian';
 export * from './shapes/2d/gaussian2D/Gaussian2D';
@@ -12,6 +13,7 @@ export type {
   PseudoVoigtShape1D,
   GaussianShape1D,
   LorentzianShape1D,
+  LorentzianDispersiveShape1D,
   GeneralizedLorentzianShape1D,
 } from './shapes/1d/Shape1D';
 export type { Shape2D, GaussianShape2D } from './shapes/2d/Shape2D';

src/shapes/1d/Shape1D.ts

@@ -18,6 +18,10 @@ export interface PseudoVoigtShape1D extends PseudoVoigtClassOptions {
   kind: 'pseudoVoigt';
 }
 
+export interface LorentzianDispersiveShape1D extends LorentzianClassOptions {
+  kind: 'lorentzianDispersive';
+}
+
 export interface GeneralizedLorentzianShape1D
   extends GeneralizedLorentzianClassOptions {
   kind: 'generalizedLorentzian';
@@ -27,4 +31,5 @@ export type Shape1D =
   | GaussianShape1D
   | LorentzianShape1D
   | PseudoVoigtShape1D
+  | LorentzianDispersiveShape1D
   | GeneralizedLorentzianShape1D;

src/shapes/1d/Shape1DInstance.ts

@@ -1,10 +1,12 @@
 import { Gaussian } from './gaussian/Gaussian';
 import { GeneralizedLorentzian } from './generalizedLorentzian/GeneralizedLorentzian';
 import { Lorentzian } from './lorentzian/Lorentzian';
+import { LorentzianDispersive } from './lorentzianDispersive/LorentzianDispersive';
 import { PseudoVoigt } from './pseudoVoigt/PseudoVoigt';
 
 export type Shape1DInstance =
   | Gaussian
   | Lorentzian
   | PseudoVoigt
+  | LorentzianDispersive
   | GeneralizedLorentzian;

src/shapes/1d/generalizedLorentzian/__tests__/GeneralizedLorentzian.shape.test.ts

@@ -1,4 +1,13 @@
-import { GeneralizedLorentzian } from '../GeneralizedLorentzian';
+import {
+  GeneralizedLorentzian,
+  calculateGeneralizedLorentzianHeight,
+  getGeneralizedLorentzianArea,
+  generalizedLorentzianFct,
+  generalizedLorentzianFwhmToWidth,
+  generalizedLorentzianWidthToFWHM,
+  getGeneralizedLorentzianFactor,
+  getGeneralizedLorentzianData,
+} from '../GeneralizedLorentzian';
 
 describe('lorentzian', () => {
   it('default factor area', () => {
@@ -20,3 +29,101 @@ describe('lorentzian', () => {
     expect(computedArea).toBeCloseTo(expectedArea, 2);
   });
 });
+
+describe('GeneralizedLorentzian class and utilities', () => {
+  it('constructor sets defaults and custom values', () => {
+    const def = new GeneralizedLorentzian();
+    expect(def.fwhm).toBe(500);
+    expect(def.gamma).toBe(0.5);
+    const custom = new GeneralizedLorentzian({ fwhm: 10, gamma: 1.5 });
+    expect(custom.fwhm).toBe(10);
+    expect(custom.gamma).toBe(1.5);
+  });
+
+  it('fwhmToWidth and widthToFWHM are inverses', () => {
+    const gLorentzian = new GeneralizedLorentzian({ fwhm: 9 });
+    const width = gLorentzian.fwhmToWidth();
+    expect(gLorentzian.widthToFWHM(width)).toBeCloseTo(9);
+  });
+
+  it('fct returns expected value for known input', () => {
+    const gLorentzian = new GeneralizedLorentzian({ fwhm: 2, gamma: 1 });
+    const val = gLorentzian.fct(0);
+    expect(typeof val).toBe('number');
+    expect(val).toBeGreaterThan(0);
+  });
+
+  it('getArea and calculateHeight are consistent', () => {
+    const gLorentzian = new GeneralizedLorentzian({ fwhm: 3, gamma: 0.5 });
+    const area = 1;
+    expect(gLorentzian.getArea(gLorentzian.calculateHeight(area))).toBe(area);
+    const height = 1;
+    expect(gLorentzian.calculateHeight(gLorentzian.getArea(height))).toBe(
+      height,
+    );
+    // area and height should be proportional
+    expect(gLorentzian.getArea(2) / gLorentzian.getArea(1)).toBe(2);
+  });
+
+  it('getGeneralizedLorentzianFactor throws for area >= 1', () => {
+    expect(() => getGeneralizedLorentzianFactor(1)).toThrow(
+      'area should be (0 - 1)',
+    );
+    expect(() => getGeneralizedLorentzianFactor(1.1)).toThrow(
+      'area should be (0 - 1)',
+    );
+  });
+
+  it('getGeneralizedLorentzianFactor returns number for area < 1', () => {
+    expect(typeof getGeneralizedLorentzianFactor(0.5)).toBe('number');
+    const gLorentzian = new GeneralizedLorentzian({ fwhm: 3, gamma: 0.5 });
+    expect(typeof gLorentzian.getFactor(0.5)).toBe('number');
+  });
+
+  it('getGeneralizedLorentzianData returns symmetric data', () => {
+    const data = getGeneralizedLorentzianData(
+      { fwhm: 5, gamma: 1 },
+      { length: 11 },
+    );
+    expect(data).toHaveLength(11);
+    for (let i = 0; i < data.length / 2; i++) {
+      expect(data[i]).toBeCloseTo(data[data.length - 1 - i]);
+    }
+  });
+
+  it('getGeneralizedLorentzianData computes length if not provided', () => {
+    const data = getGeneralizedLorentzianData({ fwhm: 2, gamma: 0.5 }, {});
+    expect(data.length % 2).toBe(1); // should be odd
+    expect(data.length).toBeGreaterThan(0);
+  });
+
+  it('getParameters returns correct parameter names', () => {
+    const gLorentzian = new GeneralizedLorentzian();
+    expect(gLorentzian.getParameters()).toStrictEqual(['fwhm', 'gamma']);
+  });
+
+  it('utility functions: fwhm/width conversion', () => {
+    const fwhm = 7;
+    const width = generalizedLorentzianFwhmToWidth(fwhm);
+    expect(generalizedLorentzianWidthToFWHM(width)).toBeCloseTo(fwhm);
+  });
+
+  it('generalizedLorentzianFct returns finite number', () => {
+    expect(Number.isFinite(generalizedLorentzianFct(0, 1, 1))).toBe(true);
+    expect(Number.isFinite(generalizedLorentzianFct(1, 2, 0.5))).toBe(true);
+  });
+
+  it('getGeneralizedLorentzianArea matches manual calculation', () => {
+    const area = getGeneralizedLorentzianArea({ fwhm: 2, height: 3, gamma: 1 });
+    expect(area).toBeCloseTo((3 * 2 * (3.14159 - 0.420894 * 1)) / 2);
+  });
+
+  it('calculateGeneralizedLorentzianHeight matches manual calculation', () => {
+    const h = calculateGeneralizedLorentzianHeight({
+      fwhm: 2,
+      area: 3,
+      gamma: 1,
+    });
+    expect(h).toBeCloseTo((3 / 2 / (3.14159 - 0.420894 * 1)) * 2);
+  });
+});

src/shapes/1d/getShape1D.ts

@@ -3,6 +3,7 @@ import { Shape1DInstance } from './Shape1DInstance';
 import { Gaussian } from './gaussian/Gaussian';
 import { GeneralizedLorentzian } from './generalizedLorentzian/GeneralizedLorentzian';
 import { Lorentzian } from './lorentzian/Lorentzian';
+import { LorentzianDispersive } from './lorentzianDispersive/LorentzianDispersive';
 import { PseudoVoigt } from './pseudoVoigt/PseudoVoigt';
 
 /**
@@ -18,6 +19,8 @@ export function getShape1D(shape: Shape1D): Shape1DInstance {
       return new Lorentzian(shape);
     case 'pseudoVoigt':
       return new PseudoVoigt(shape);
+    case 'lorentzianDispersive':
+      return new LorentzianDispersive(shape);
     case 'generalizedLorentzian':
       return new GeneralizedLorentzian(shape);
     default: {

src/shapes/1d/lorentzian/Lorentzian.ts

@@ -10,7 +10,7 @@ export interface LorentzianClassOptions {
   fwhm?: number;
 }
 
-interface GetLorentzianAreaOptions {
+export interface GetLorentzianAreaOptions {
   /**
    * The maximum intensity value of the shape
    * @default 1

src/shapes/1d/lorentzianDispersive/LorentzianDispersive.ts

@@ -0,0 +1,86 @@
+import type { GetData1DOptions } from '../GetData1DOptions';
+import type { Parameter, Shape1DClass } from '../Shape1DClass';
+import {
+  calculateLorentzianHeight,
+  getLorentzianFactor,
+  LorentzianClassOptions,
+  lorentzianFwhmToWidth,
+  lorentzianWidthToFWHM,
+} from '../lorentzian/Lorentzian';
+
+export class LorentzianDispersive implements Shape1DClass {
+  /**
+   * Full width at half maximum.
+   * @default 500
+   */
+  public fwhm: number;
+
+  public constructor(options: LorentzianClassOptions = {}) {
+    const { fwhm = 500 } = options;
+
+    this.fwhm = fwhm;
+  }
+
+  public fwhmToWidth(fwhm = this.fwhm) {
+    return lorentzianFwhmToWidth(fwhm);
+  }
+
+  public widthToFWHM(width: number) {
+    return lorentzianWidthToFWHM(width);
+  }
+
+  public fct(x: number) {
+    return lorentzianDispersiveFct(x, this.fwhm);
+  }
+
+  //eslint-disable-next-line
+  public getArea(_height: number) {
+    return 0;
+  }
+
+  public getFactor(area?: number) {
+    return getLorentzianFactor(area);
+  }
+
+  public getData(options: GetData1DOptions = {}) {
+    return getLorentzianDispersiveData(this, options);
+  }
+
+  public calculateHeight(area = 1) {
+    return calculateLorentzianHeight({ fwhm: this.fwhm, area });
+  }
+
+  public getParameters(): Parameter[] {
+    return ['fwhm'];
+  }
+}
+
+export const lorentzianDispersiveFct = (x: number, fwhm: number) => {
+  return (2 * fwhm * x) / (4 * x ** 2 + fwhm ** 2);
+};
+
+export const getLorentzianDispersiveData = (
+  shape: LorentzianClassOptions = {},
+  options: GetData1DOptions = {},
+) => {
+  let { fwhm = 500 } = shape;
+  let {
+    length,
+    factor = getLorentzianFactor(),
+    height = calculateLorentzianHeight({ fwhm, area: 1 }),
+  } = options;
+
+  if (!length) {
+    length = Math.min(Math.ceil(fwhm * factor), Math.pow(2, 25) - 1);
+    if (length % 2 === 0) length++;
+  }
+
+  const center = (length - 1) / 2;
+  const data = new Float64Array(length);
+  for (let i = 0; i <= center; i++) {
+    data[i] = lorentzianDispersiveFct(i - center, fwhm) * height;
+    data[length - 1 - i] = -data[i];
+  }
+
+  return data;
+};

src/shapes/1d/lorentzianDispersive/__tests__/LorentzianDispersive.fct.test.ts

@@ -0,0 +1,10 @@
+import { LorentzianDispersive } from '../LorentzianDispersive';
+
+describe('Lorentzian function', () => {
+  it('Lorentzian.fct', () => {
+    const lorentzian = new LorentzianDispersive({ fwhm: 0.2 });
+    expect(lorentzian.fct(0)).toBeCloseTo(0);
+    expect(lorentzian.fct(0.1)).toBeCloseTo(0.5);
+    expect(lorentzian.fct(-0.1)).toBeCloseTo(-0.5);
+  });
+});

src/shapes/1d/lorentzianDispersive/__tests__/LorentzianDispersive.shape.test.ts

@@ -0,0 +1,76 @@
+import { ROOT_THREE } from '../../../../util/constants';
+import { LorentzianDispersive } from '../LorentzianDispersive';
+
+describe('lorentzian', () => {
+  it('default factor area', () => {
+    const lorentzian = new LorentzianDispersive();
+    const data = lorentzian.getData();
+    expect(data).toHaveLength(3183099);
+    const area = data.reduce((a, b) => a + b, 0);
+    expect(area).toBeCloseTo(0, 5);
+    const expectedArea = 0;
+    const height = lorentzian.calculateHeight(expectedArea);
+    const computedArea = lorentzian.getArea(height);
+    expect(computedArea).toBeCloseTo(expectedArea, 2);
+    expect(lorentzian.getParameters()).toStrictEqual(['fwhm']);
+
+    const width = 20;
+    expect(lorentzian.widthToFWHM(width)).toBe(width * ROOT_THREE);
+    expect(lorentzian.fwhmToWidth(lorentzian.widthToFWHM(width))).toBe(width);
+  });
+  it('differents factor', () => {
+    const lorentzian = new LorentzianDispersive({ fwhm: 1000 });
+    //areas coverage by the real part of the lorentzian
+    const areas = [0.98, 0.96, 0.7, 0.4, 0.2];
+    for (let area of areas) {
+      const data = lorentzian.getData({ factor: lorentzian.getFactor(area) });
+      const sum = data.reduce((a, b) => a + b, 0);
+      expect(sum).toBeCloseTo(0, 3);
+    }
+  });
+  it('default factor', () => {
+    const fwhm = 10;
+    const lorentzian = new LorentzianDispersive({ fwhm });
+    const data = lorentzian.getData({ height: 1 });
+    expect(data).toHaveLength(63663);
+    const center = (data.length - 1) / 2;
+    expect(Math.abs(data[center])).toBe(0);
+    expect(data[center + fwhm / 2]).toBe(0.5);
+    expect(data[center - fwhm / 2]).toBe(-0.5);
+    const area = data.reduce((a, b) => a + b, 0);
+    expect(area).toBeCloseTo(0, 5);
+  });
+  it('fwhm 10, factor 500', () => {
+    const fwhm = 10;
+    const lorentzian = new LorentzianDispersive({ fwhm });
+    const length = 5001;
+    const data = lorentzian.getData({ length });
+    const height = lorentzian.calculateHeight();
+    const center = Math.floor((length - 1) / 2);
+    expect(Math.abs(data[center])).toBe(0);
+    expect(data[center - fwhm / 2]).toBe(-height / 2);
+    expect(data[center + fwhm / 2]).toBe(height / 2);
+    expect(data).toHaveLength(length);
+    const area = data.reduce((a, b) => a + b, 0);
+    expect(area).toBeCloseTo(0, 2);
+  });
+  it('odd fwhm', () => {
+    const lorentzian = new LorentzianDispersive({ fwhm: 11 });
+    const data = lorentzian.getData({ length: 11, height: 2 });
+    const lenG = data.length;
+    const center = Math.floor((lenG - 1) / 2);
+    expect(data[center]).toBeCloseTo(0, 4);
+    expect(data[center - 1]).toBeCloseTo(-data[center + 1], 4);
+    expect(data[center]).toBeGreaterThan(data[center - 1]);
+    expect(data[center]).toBeLessThan(data[center + 1]);
+  });
+  it('even fwhm', () => {
+    const lorentzian = new LorentzianDispersive({ fwhm: 10 });
+    const data = lorentzian.getData({ length: 10, height: 1 });
+    const lenG = data.length;
+    const center = Math.floor((lenG - 1) / 2);
+    expect(data[center]).toBeCloseTo(-data[center + 1], 4);
+    // in the infinity both should be zero but in the practice:
+    expect(data[0]).toBeCloseTo(-data[data.length - 1], 4);
+  });
+});