Skip to content

[bug] PROFILE_LORENTZ uses the wrong sign for Delta0 #57

@cexen

Description

@cexen

In PROFILE_LORENTZ, the pressure shift Delta0 has the wrong sign.

Current implementation:

def PROFILE_LORENTZ(Nu,Gamma0,Delta0,WnGrid,YRosen=0.0,Sw=1.0):
    if YRosen==0.0:
        return Sw*Gamma0/(pi*(Gamma0**2+(WnGrid+Delta0-Nu)**2))
    else:
        return Sw*(Gamma0+YRosen*(WnGrid+Delta0-Nu))/(pi*(Gamma0**2+(WnGrid+Delta0-Nu)**2))

Expected (sign of Delta0 flipped):

def PROFILE_LORENTZ(Nu,Gamma0,Delta0,WnGrid,YRosen=0.0,Sw=1.0):
    if YRosen==0.0:
        return Sw*Gamma0/(pi*(Gamma0**2+(WnGrid-Delta0-Nu)**2))
    else:
        return Sw*(Gamma0+YRosen*(WnGrid-Delta0-Nu))/(pi*(Gamma0**2+(WnGrid-Delta0-Nu)**2))

Minimal Reproduction:

import hapi
import numpy as np

TABLE = "H₂¹⁶O"
hapi.db_begin("hapi_db")
hapi.fetch(TABLE, 1, 1, 3744.65, 3744.66)

nu = hapi.getColumn(TABLE, "nu")
delta = hapi.getColumn(TABLE, "delta_air")
assert len(nu) == 1
nu0 = nu[0]
delta_air = delta[0]
print(f"nu0 = {nu0}")
print(f"delta_air = {delta_air}")


def _run(fn, label):
    nu_out, coef = fn(
        SourceTables=TABLE,
        WavenumberRange=(nu0 - 0.05, nu0 + 0.05),
        WavenumberStep=0.000001,
    )
    peak = nu_out[np.argmax(coef)]
    shift = peak - nu0
    print(f"{label}: {peak:.6f} ({shift:+.6f})")


_run(hapi.absorptionCoefficient_Priority, "Priority")
_run(hapi.absorptionCoefficient_HT, "HT")
_run(hapi.absorptionCoefficient_SDVoigt, "SDVoigt")
_run(hapi.absorptionCoefficient_Voigt, "Voigt")
_run(hapi.absorptionCoefficient_Lorentz, "Lorentz")
_run(hapi.absorptionCoefficient_Doppler, "Doppler")

Output (excerpt):

nu0 = 3744.651222
delta_air = -0.005164

Priority: 3744.646058 (-0.005164)
HT: 3744.646058 (-0.005164)
SDVoigt: 3744.646058 (-0.005164)
Voigt: 3744.646058 (-0.005164)
Lorentz: 3744.656386 (+0.005164)
Doppler: 3744.651222 (-0.000000)

Lorentz shifts in the opposite direction.

Metadata

Metadata

Assignees

No one assigned

    Labels

    No labels
    No labels

    Projects

    No projects

    Milestone

    No milestone

    Relationships

    None yet

    Development

    No branches or pull requests

    Issue actions