Performing Fermi Edge Corrections

The most straightforward way to correct the Fermi edge, either due to monochromator miscalibration in a photon energy scan or due to using a straight slit on a hemispherical analyzer is just to broadcast an edge.

In the case of correcting for the slit shape, it may be helpful to further fit a model for the edge shape, like a quadratic, so that a smooth correction is applied across the detector.

Correcting monochromator miscalibration

from lmfit.models import ConstantModel
import arpes
from arpes.fits import AffineBroadenedFD
from arpes.io import example_data

photon_energy = example_data.photon_energy

Activating auto-logging. Current session state plus future input saved.
Filename       : logs/unnamed_2026-03-24_23-29-47.log
Mode           : backup
Output logging : False
Raw input log  : False
Timestamping   : False
State          : active

edge_data = photon_energy.sel(phi=slice(-0.28, -0.15), eV=slice(-0.1, 0.1)).sum("phi").spectrum
edge_data.T.plot()

<matplotlib.collections.QuadMesh at 0x7cb931f03320>

import matplotlib.pyplot as plt

model = AffineBroadenedFD() + ConstantModel()
# fermi_param = AffineBroadenedFD().guess(data=edge_data.sel(hv=70, method="nearest").values, x=edge_data.coords["eV"].values)

param = AffineBroadenedFD().make_params(
    center=0.0, width=0.02, sigma=0.02, lin_slope=0, const_bkg=15000
)
fit_results = edge_data.S.modelfit("eV", model, params=param)
# fit_results = broadcast_model(AffineBroadenedFD, edge_data, "hv")

fig, ax = plt.subplots()
edge_data.T.plot(ax=ax)
# ax.scatter(*fit_results.results.F.p("center").G.to_arrays(), color="red")
ax.scatter(
    fit_results.modelfit_results.F.p("center").coords["hv"],
    fit_results.modelfit_results.F.p("center").values,
    color="red",
)

<matplotlib.collections.PathCollection at 0x7cb931dd2000>

Now, we can perform the shift to correct the edge.

corrected_photon_energy = photon_energy.spectrum.G.shift_by(
    fit_results.modelfit_results.F.p("center"), shift_axis="eV", shift_coords=True
)  # Note that G.shift_by is applied to "spectrum" (i.e. xr.DataArray)

corrected_edge_data = corrected_photon_energy.sel(phi=slice(-0.28, -0.15),
                                                  eV=slice(-0.1, 0.1)).sum(
    "phi"
)

model = AffineBroadenedFD() + ConstantModel()
param = AffineBroadenedFD().make_params(
    center=-0.025, width=0.02, sigma=0.02, lin_slope=0, const_bkg=15000,
)
results_check = corrected_edge_data.S.modelfit("eV", model, params=param)

fig, ax = plt.subplots()
corrected_edge_data.T.plot(ax=ax)
# ax.scatter(*results_check.results.F.p("center").G.to_arrays(), color="red")
ax.scatter(
    results_check.modelfit_results.F.p("center").coords["hv"],
    results_check.modelfit_results.F.p("center").values,
    color="red",
)

<matplotlib.collections.PathCollection at 0x7cb931448f50>

As we can see, the edge is now uniform across photon energy and correctly zero referenced.

A caveat to be aware of when shifting data is whether to make data adjustments only or to use coordinate adjustments as well. Coordinate adjustments (above: shift_coords=True) are useful when the shift is very large. If the coords are not allowed to compensate for some of the shift in that context, large portions of data will be shifted out of the array extent and be replaced by np.nan or 0.

However, if coordinates no longer agree between two pieces of data, we will not be able to perform array operations involving both of them, because of their incompatible coordinates.

The correct behavior is context dependent and requires you to consider what analysis you are trying to do.

Correcting the curved Fermi edges

Let’s now turn to an example addressing the uneven energy calibration arising from the use of a straight slit in ARPES data.

We will shortly turn to the question of momentum conversion, but we will want to have this issue corrected before converting. The correction we need to apply is a function of the detector angle phi, so it will not be a constant function of any momentum coordinate, in general.

First, let’s assess the data.

from arpes.io import example_data

cut = example_data.map.sum("theta").spectrum
cut = cut.sel(eV=slice(-0.2, 0.1), phi=slice(-0.25, 0.3))
cut.S.plot()

Now, let’s fit edges to this data.

import matplotlib.pyplot as plt
from lmfit.models import QuadraticModel, ConstantModel
from arpes.fits.fit_models import AffineBroadenedFD

model = AffineBroadenedFD() + ConstantModel()
params = AffineBroadenedFD().make_params(
    center=0,
    width=0.005,
    sigma=0.02,
    const_bkg=200000,
    lin_slope=0,
)


fit_results = cut.S.modelfit("eV", model, params=params)

fit_results.modelfit_results.F.plot_param("center")
plt.gca().set_ylim([-0.05, 0.05])

(-0.05, 0.05)

We could either refine these fits a little by setting some constraints, or we can make a smooth correction by fitting a quadratic to these edge locations.

quad_mod = (
    fit_results.modelfit_results.F.p("center")
    .S.modelfit("phi", QuadraticModel())
    .modelfit_results.item()
)

# quad_mod = QuadraticModel().guess_fit(fit_results.results.F.p("center"))
quad_mod_plot = quad_mod.plot()

Now, we can shift the data exactly like we did before. In order to get the shift amount to apply at each phi, we just evaluate our quadratic at these phi.

fmap = example_data.map.spectrum
edge = quad_mod.eval(x=fmap.phi)

corrected_map = fmap.G.shift_by(edge, shift_axis="eV", by_axis="phi")
extend_corrected_map = fmap.G.shift_by(edge,
                                       shift_axis="eV",
                                       by_axis="phi",
                                       extend_coords=True)

fig, axes = plt.subplots(1, 3, figsize=(10, 5))
corrected_map.isel(theta=10).S.plot(ax=axes[0])
fmap.isel(theta=10).S.plot(ax=axes[1])
extend_corrected_map.isel(theta=10).S.plot(ax=axes[2])
axes[0].set_title("Corrected data")
axes[1].set_title("Original data")
axes[2].set_title("Corrected data, coordinate extended")
fig.tight_layout()

Much better. With this in order, we can now consider momentum conversion.

Exercises

1. Correct the Fermi edge on the Bi2Se3 cut data example_data.cut. How does the presence of the dispersing surface state affect the range you should use?

2. Perform a random shift of a Fermi edge using .G.shift_by. Then, try to correct it. How close is your recovered correction to the original shift?