No, that has to be a different definition of "exponential." The plot in this post
(https://smf.probesoftware.com/index.php?msg=11246) shows the PfE-implemented exponential deadtime correction is larger than the logarithmic correction, but my exponential correction is smaller than the logarithmic correction. Also, the PfE correction is shown to have a discontinuity at 1/e, before the logarithmic correction's discontinuity at 1 - 1/e (approx 0.631). 1/(2-exp^x) has its discontinuity at ln(2) (approx 0.693), which is further right than the logarithmic correction.

Quote from: aducharme on April 29, 2026, 11:42:12 PMReading the deadtime correction paper, the logarithm appears to have fallen out by extrapolating the x + x^2/2 terms as the first two terms in the Taylor series of log(1-x)=\sum x^n/n. To me, those first two terms are hallmarks of the Taylor series for e^x = \sum x^n/n!. This gives another possible deadtime correction i/(2-e^(i*\tau)). It is surprisingly similar to the logarithm-based correction:
Maybe this would get rid of the squiggle? Or perhaps it's not a good model. John, would you be able to give it a shot?

Yes, we implemented this dead time correction, but the problem is that it "blows up" at relatively low count rates depending on the magnitude of the dead time constant. See the last dead time option here in Probe for EPMA:

Quote from: JonF on April 30, 2026, 01:59:03 AM

Quote from: Ben Buse on April 29, 2026, 09:11:05 AMJohn, what do you think about loading in specified concentration maps into CalcImage? If the EDS maps were quantified elsewhere, and we loaded the maps in as specified elements to allow for correct matrix correction of the WDS elements?

Just to add my two pennies worth, the problem with this is that you would need to quant the EDS maps first, meaning that for these to be accurate, you need to quantify all of the major elements. This would be fine if you were just using the WDS to look at minor/trace phases (alongside the EDS), but if you're using the WDS to disentangle a nasty interference of major elements, then you'd potentially be loading in the wrong quant EDS data, whether by omission of the WDS-measured elements (and hence incorrect matrix corrections), or by including elements that required the WDS in the first place (and hence have a major interference issue).

Jon is exactly correct.

In addition to *all* elements needing to be included in the matrix corrections (and dynamic specified elements mentioned above), if spectral interferences are present, we also need to have both the interfered *and* interfering elements elements present for the quantitative interference corrections.

These compositionally dependent corrections also includes other corrections such as the mean atomic number (background) correction and excess oxygen calculations from ferric iron (which requires all cations be measured/specified).

Here's a modification on the load script of Yuji to trim the first 160 channels of the spectrum, so that the load spectrum matches the emsa spectra of the whole map




import os
import numpy as np
#import DigitalMicrograph as DM

#dir=r"C:\Users\arksa\Documents\em-data\jeol-pts"
#dir=r"C:\Users\glxbb\OneDrive - University of Bristol\JEOL8530F\EpmaData\Ben\eds wds sjio pfe\eds wds sjio pfe_0004_MAP\Pos_0002\EDS"
#dir=r"C:\Users\glxbb\OneDrive - University of Bristol\JEOL8530F\EpmaData\Ben\eds test2\eds test2_0001_MAP\Pos_0003\EDS"
dir=r"C:\Users\glxbb\OneDrive - University of Bristol\JEOL8530F\EpmaData\Ben\eds wds map ga garnet\eds wds map ga garnet_0001_MAP\Pos_0001\EDS"
os.chdir(dir)
file=r"2406PTTD.PTS"

x=60
y=60

# number of energy channel of EDS.
ch=4096

# binning of energy channel in output.
# if you want enegy resolution, set 1.
# this parameter is set just for memory and disc saving.
binn=1

# cutoff energy of output.
# this is also for memory saving.
cut=2048-160

data=np.fromfile(file,dtype="u2",offset=16384)

flag=np.arange(data.size)[(data==55296) | (data==53504) | (data==53760)]
# assuming livetime stamp is 24576 = 0x6000
# assuming realtime stamp is 28672 = 0x7000
livet=24576
realt=28672

# resolution of time stamp.
# maybe 10 msec
tscale=10

#data=np.where((data==24576) | (data==28672), 1, data-45056)
data=np.where((data==livet) | (data==realt), binn*((data-livet)//(realt-livet)), data-45056)


flagsize=flag.size

print(flagsize)

si=np.zeros((x+1,y,(ch-cut-160)//binn),dtype="u2")
siTest=np.zeros((x+1,y,(ch-cut)//binn),dtype="u2")

ypos=-1
j=0
for i in range(flag.size):
#for i in range(10):
    pos=flag[i]
    if data[pos]==55296-45056:
        ypos+=1
        j=0
    elif data[pos]==53760-45056:
        nextpos=flag[i+1]
        if data[nextpos]==53760-45056 or data[nextpos]==53504-45056:
            si[j,ypos,:]=np.bincount(data[pos+1:nextpos]//binn,minlength=(ch-cut)//binn)[:(ch-cut)//binn][160:]
            #si[j,ypos,:]=np.bincount(data[pos+1:nextpos])[160:2048+160]
            #SiStore=np.bincount(data[pos+1:nextpos])
            j+=1

si[:,:,0:2]*=tscale

si2=np.transpose(si,(2,0,1))
#siDM=DM.CreateImage(si2.copy(order="C"))
#siDM.ShowImage()
#del siDM




There is a recalibrate function that might work to change the energy offset.
"""
    recalibrate(s::Spectrum{T}, es::LinearEnergyScale)

Allows changing the energy scale on a spectrum from one LinearEnergyScale to another as though the spectrum were
measured on a different detector.  The algorithm uses a FFT-base scheme to rescale and shift the spectral data.
This scheme allows for fractional shifts of offset and fractional changes in the width.  It is limited in that
the change in width must produce an integral number of channels in the resulting spectrum.  The algorithm
maintains the total spectrum integral so the new spectrum can be used for quantitative purposes.
Plotting one spectrum over the other should maintain peak position but is likely to change the channel counts.
"""
function recalibrate(s::Spectrum{T}, es::LinearEnergyScale)

Example:
les = LinearEnergyScale(3.0, 10.1)
spec = loadspectrum("\path\to\spectrum\spec.msa")
new_spec = recalibrate(spec, les)


Quote from: Nicholas Ritchie on April 30, 2026, 06:00:34 AM

Quote from: Ben Buse on April 30, 2026, 03:48:51 AMCan I make a map sum a standard?

Code Select Expand

MapSum = sum(hs3)
MapSum[:LiveTime]=hs3[1,1][:LiveTime]*60*60
maprefs = references([reference(n"Al",MapSum,mat"Ca3Mn57Al40Si60O240")],128)
ERROR: DomainError with -2082.385000000002:
sqrt was called with a negative real argument but will only return a complex result if called with a complex argument. Try sqrt(Complex(x)).


I'm not exactly sure where the problem lies.  It may be a overflow in

Code Select Expand

hs3[1,1][:LiveTime]*60*60
Try:

Code Select Expand

float(hs3[1,1][:LiveTime])*60*60  # or
hs3[1,1][:LiveTime]*60.0*60.0


Here the problem seems to be using the map spectra as a reference

maprefs = references([reference(n"Al",hs3[1,1],mat"Al")],128)
ERROR: DomainError with -2082.385000000002:
sqrt was called with a negative real argument but will only return a complex result if called with a complex argument. Try sqrt(Complex(x)).
Stacktrace:
  [1] throw_complex_domainerror(f::Symbol, x::Float64)
    @ Base.Math .\math.jl:33
  [2] sqrt
    @ .\math.jl:627 [inlined]
  [3] resolution_to_fwhm
    @ C:\Users\glxbb\.julia\packages\NeXLSpectrum\eVZfj\src\detector.jl:155 [inlined]
  [4] resolution
    @ C:\Users\glxbb\.julia\packages\NeXLSpectrum\eVZfj\src\detector.jl:165 [inlined]
  [5] resolution(eV::Float64, det::BasicEDS)
    @ NeXLSpectrum C:\Users\glxbb\.julia\packages\NeXLSpectrum\eVZfj\src\detector.jl:309
  [6] buildfilter(::Type{Float64}, ty::Type{G2Filter}, det::BasicEDS, a::Float64, b::Float64)
    @ NeXLSpectrum C:\Users\glxbb\.julia\packages\NeXLSpectrum\eVZfj\src\filter.jl:320
  [7] buildfilter
    @ C:\Users\glxbb\.julia\packages\NeXLSpectrum\eVZfj\src\filter.jl:312 [inlined]
  [8] references(refs::Vector{NeXLSpectrum.ReferencePacket}, det::BasicEDS; ftype::Type{Float64}, filter::Type)
    @ NeXLSpectrum C:\Users\glxbb\.julia\packages\NeXLSpectrum\eVZfj\src\reference.jl:156
  [9] references
    @ C:\Users\glxbb\.julia\packages\NeXLSpectrum\eVZfj\src\reference.jl:146 [inlined]
 [10] references(refs::Vector{NeXLSpectrum.ReferencePacket}, fwhm::Int64; kwargs::@Kwargs{})
    @ NeXLSpectrum C:\Users\glxbb\.julia\packages\NeXLSpectrum\eVZfj\src\reference.jl:166
 [11] top-level scope
    @ REPL[255]:1

As I post just above I'm wondering how to use the emsa standards when map has energy offset

Thanks Nicholas, I will give it a go.

Now I'm trying to figure out how to make my map and standards have the same detector, when the map has a large energy offset, compared to the emsa standards which have no energy offset. I've reloaded map with 2048 channels

60 × 60 HyperSpectrum{UInt16,(:Y, :X)}[8 sp1 spess, -1638.3 + 10.0⋅ch eV, 2048 ch]
refs4
References[
        BasicEDS[2048 chs, 0.0 + 10.0⋅ch eV, 130.0 eV @ Mn K-L3, 2 ch LLD, [Be,Sc,Ba,Pu]]

so
resmap4 = fit_spectrum(hs4,refs4)
ERROR: AssertionError: The detector for the hyper-spectrum must match the detector for the filtered references.
Stacktrace:

Any suggestions

Quote from: Ben Buse on April 30, 2026, 03:48:51 AMCan I make a map sum a standard?

Code Select Expand

MapSum = sum(hs3)
MapSum[:LiveTime]=hs3[1,1][:LiveTime]*60*60
maprefs = references([reference(n"Al",MapSum,mat"Ca3Mn57Al40Si60O240")],128)
ERROR: DomainError with -2082.385000000002:
sqrt was called with a negative real argument but will only return a complex result if called with a complex argument. Try sqrt(Complex(x)).


I'm not exactly sure where the problem lies.  It may be a overflow in
hs3[1,1][:LiveTime]*60*60
Try:
float(hs3[1,1][:LiveTime])*60*60  # or
hs3[1,1][:LiveTime]*60.0*60.0


You can also construct references in two steps.  1) Load the spectrum and set properties. 2) Replace the reference spectrum path with the loaded spectrum.

spec=pathspectrum("path\to\references\SiO2")
spec[:LiveTime] = 60.0
reference(n"Si", spec, mat"SiO2")


Quote from: Ben Buse on April 30, 2026, 01:11:21 AMLooking at this

Code Select Expand

LiveTime => lt,Seems live time is fixed for every pixel in the map is it possible to use an array of live time

Using Nicholas example above, I loaded the map save from hyperspy as rpl,

Code Select Expand

hs = NeXLSpectrum.compress(HyperSpectrum(
                  LinearEnergyScale(0.0,40.0),
                  Dict{Symbol,Any}(
                    :TakeOffAngle => deg2rad(35.0),
                                :ProbeCurrent => 1.0,
                    :LiveTime => lt,
                    :BeamEnergy => 30.0e3,
                    :Name => "8 sp1 spess"
                  ),
                  readrplraw(raw"C:\\Users\\glxbb\\julia\\edstry2"),
                            fov = [ 60u"mm", 60u"mm"], offset= [ 0.0u"mm", 0.0u"mm" ]
              ))
60 × 60 HyperSpectrum{UInt16,(:Y, :X)}[8 sp1 spess, 0.0 + 40.0⋅ch eV, 512 ch]


You can specify per pixel live-times by adding the `livetime` argument to the HyperSpectrum constructor assigning it a matrix the same dimensions as the RPL/RAW data. (Remove the `:LiveTime => lt,`)
    readrplraw(....),
    livetime=fill(get(props, :LiveTime, 1.0), (512, 512))
)