Spectrum

Microanalytical X-ray Spectrum Analysis

Spectrum Manipulation

See the Spectrum Methods page for the most used methods and details.

NeXLSpectrum.duane_huntFunction
duane_hunt(spec::Spectrum)

Estimates the Duane-Hunt limit (the energy at which the continuum goes to effectively zero.)

source
NeXLSpectrum.sigmaFunction
sigma(spec::Spectrum, specs::AbstractArray{<:Spectrum}, chs::AbstractRange{<:Integer})::Vector{Float64}

Computes on a channel-by-channel basis how much spec spectrum deviates from the mean of the other spectra in specs. The result is expressed in terms of the standard deviation expected from count statistics alone. Assuming spec varies only by count statistics we expect the result values have a mean 0.0 and a standard deviation of 1.0.

source
NeXLSpectrum.findsimilarFunction
findsimilar(
    specs::AbstractArray{Spectrum{T}}; 
    atol = nothing, 
    rtol=1.5, 
    minspecs=3
)::Vector{Spectrum{T}}
findsimilar(
    specs::AbstractArray{Spectrum{T}},
    det::Detector,
    elm::Element; 
    atol = nothing, 
    rtol=1.5, 
    minspecs = 3
)::Vector{Spectrum{T}}

Filters a collection of spectra for the ones most similar to the average by removing the least similar spectrum sequentially until all the remaining spectra are either:

  • less than atol (if atol != nothing)
  • less than rtol * median(others) (if rtol != nothing)

when applying the 'similarity(...)` function to the spectrum and the sum of the other spectra.

This is useful for finding which of a set of replicate spectra are sufficiently similar to each other.

source
NeXLSpectrum.multiscoreFunction
multiscore(specs::AbstractArray{<:Spectrum}, e0=specs[1][:BeamEnergy])

Compares each spectrum against the mean spectrum by comparing the 200 eV to 500 eV range with the [e0/2,e0/1.5] range. If all spectra are equivalent at low energies then all the scores will be close to zero. A number less than zero means that the low energy region of the spectrum was depressed relative to the others. A number more than zero means that the high energy region of the spectrum was elevated relative to the others. The multi-score is sensitive to tilt and obstructions like surface texture which may make one spectrum's low energy be more absorbed than the others.

source
NeXLSpectrum.multirankFunction
multirank(specs::AbstractArray{<:Spectrum})::Float64

A single number that compares the low and high energy portions of the spectra for similarity. A multirank(...) score of zero means all the spectra are very similar and a large number means very different. A high score suggests that one or more of the spectra may suffer from additional low energy absorption due to surface roughness, an obstruction, sample tilt or other.

source
NeXLSpectrum.plot_compareFunction
plot_compare(specs::AbstractArray{<:Spectrum}, mode=:Plot; xmin=100.0, xmax=1.0, palette = NeXLPalette)

Plots a comparison of the channel-by-channel data from each individual spectrum relative to the dose-corrected mean of the other spectra. Count statistics are taken into account so if the spectra agree to within count statistics we expect a mean of 0.0 and a standard deviation of 1.0 over all channels. Note: xmax is relative to the :BeamEnergy.

source
NeXLSpectrum.applyFunction
apply(filt::SavitzkyGolayFilter, spec::Spectrum, applyLLD=false)

Applys a function to the channel data in spec (with/wo the low-level discriminator.) The function can only be a function of the counts data and can not change the energy scale or spectrum properties. The result is a Spectrum.

Example:

apply(SavitzkyGolayFilter{6,4}(),spec)
source
NeXLSpectrum.loadmultispecFunction
loadmultispec(path::AbstractString, basefn::AbstractString; indexes=0:3, fnmapper::String = "{1}[{2}].msa")

Load multiple spectra using the basefn and fnmapper to determine which spectra to load. The spectra should be related in the sense that they were all collected simulataneously so they have the same :RealTime, :BeamEnergy and :LiveTime.

source
NeXLUncertainties.uvFunction
uv(spec::Spectrum, chs::AbstractRange{<:Integer}=eachindex(spec))::Vector{UncertainValue}

Converts the count's data in a spectrum into an Vector{UncertainValue} assuming count statistics can be approximated by C ± √C.

source
NeXLSpectrum.χ²Function
χ²(s1::Spectrum{T}, s2::Spectrum{U}, chs)::T where {T<:Real, U <: Real}
χ²(specs::AbstractArray{Spectrum{T}}, chs)::Matrix{T}

Computes the dose corrected reduced χ² metric between s1 and s2 over the channels in chs.

The second form computes a matrix of χ² comparing each spectrum in the array to the others.

source
NeXLSpectrum.recalibrateFunction
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 algorith uses a FFT-base scheme to rescale and shift the spectral data. Ths 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.

source
NeXLSpectrum.shiftFunction
shift(s::Spectrum, ev::AbstractFloat)::Spectrum

Shift the entire spectrum along the energy axis by a specified number of ev by modifying the counts data. This function can shift by a fractional number of channels.

source
NeXLSpectrum.offsetFunction
offset(s::Spectrum, dcounts::Real)

Returns a Spectrum like s but with dcounts added to each channel.

source
NeXLSpectrum.dosenormalizeFunction
dosenormalize(spectrum::Spectrum{T}, dose=60.0)::Spectrum{T} where { T <: Real }
dosenormalize(spectrum::Spectrum{T}, dose=60.0)::Spectrum{Float64} where { T <: Integer }

Compute a spectrum which is spectrum rescaled to a live time times probe current equal to dose. Useful for setting spectra on an equivalent acquisition duration scale.

source
NeXLSpectrum.extentFunction
extent(xrayE::Float64, res::Resolution, ampl::Float64)

Calculates the extent of the peak interval for an x-ray of the specified energy.

source
extent(escape::EscapeArtifact, res::Resolution, ampl::Float64)::Tuple{2,Float64}

The extent of an escape artifact is determined by the resolution of the detector at the energy of the escape peak.

source
extent(escape::ComptonArtifact, res::Resolution, ampl::Float64)::Tuple{2,Float64}

The extent of a Compton artifact is determined by the resolution of the detector at the energy of the Compton peak.

source
extent(sf::SpectrumFeature, det::Detector, ampl::Float64)::Tuple{Float64, Float64}

Computes the channel range encompassed by the specified set of x-ray transitions down to an intensity of ampl. Relative line weights are taken into account.

source
NeXLSpectrum.characteristiccountsFunction
characteristiccounts(ffr::FiterFitResult, strip)

Number of spectrum counts that were accounted for by the fitted elements with the strip Element(s) removed.

source
NeXLSpectrum.sumcountsFunction
sumcounts(hss::HyperSpectrum, cis::CartesianIndices = CartesianIndices(hss))

An array containing the number of counts at each pixel.

source
NeXLSpectrum.shannon_entropyFunction
shannon_entropy(spec::Spectrum)

Computes a measure of the information content in a spectrum. As there become more and more distinct values in a spectrum, this value approaches log2(nchannels(spec)). This number reflects the number of bits necessary to encode the spectrum data with maximum efficiency.

This is inspired by John Colby's FLAME software which did something similar. Although, to be honest, I don't know how his algorithm was implemented.

source
shannon_entropy(img::AbstractArray{Gray{N0f8}})

Computes the log2-entropy of the data in the image. The entropy(...) in Images.jl 24.1 is buggy and is removed in 25.0

source
NeXLSpectrum.similarityFunction
similarity(s1::Spectrum{T}, s2::Spectrum{T}, chs)::Float64 where {T<:Real}
similarity(specs::AbstractArray{Spectrum{T}}, chs)::Vector{Float64}
similarity(specs::AbstractArray{Spectrum}, minE::Float64=100.0)::Vector{Float64}
similarity(specs::AbstractArray{<:Spectrum}, det::Detector, elm::Element)::Vector{Float64}
similarity(specs::AbstractArray{Spectrum{T}}, det::Detector, mat::Material)::Vector{Float64}

Returns a vector of similarity metrics which measure how similar the i-th Spectrum is to the other spectra. The mean reduced χ² statistic metric is such that if s1 and s2 differ by only count statistics then the metric will be approximately unity. If s1 and s2 vary due to probe current drift, sample inhomogeneity, surface roughness or other non-count statistics related reasons then the metric will be larger than one.

The first version covers all the channels between minE and the nominal beam energy. The third and fourth versions considers those channels representing peaks in a spectrum from the Material or Element on the Detector.

source
NeXLSpectrum.fittedcontinuumFunction

fittedcontinuum( spec::Spectrum, det::EDSDetector, resp::AbstractArray{<:Real,2}; # mode = :Global [ | :Local ] # Fit to all ROIs simultaneously (:Global) or to each roi independently (:Local) minE::Float64 = 1.5e3, maxE::Float64 = 0.95 * spec[:BeamEnergy], width::Int = 20, # Width of ROI at each end of each patch of continuum that is matched brem::Type{<:NeXLBremsstrahlung} = Castellano2004a, mc::Type{<:MatrixCorrection} = Riveros1993, )::Spectrum

Fit the continuum under the characteristic peaks by fitting the closest continuum ROIs. The low energy peaks are fit using the continuum immediately higher in energy and the high energy peaks are fit using the continuum on both sides.

  • mode = :Global [ | :Local ] Global fits the model to the data using a single scale factor. :Local tweaks the

global model at ROIs above and below the characteristic peaks.

Typically, :Global produces the overall best fit but :Local fits better around the characteristic peaks and is better for modeling the continuum under the peaks.

source

Spectrum Plotting

Gadfly.plotMethod
plot(
    specs::Union{Spectrum...,AbstractVector{Spectrum{<:Real}}};
    klms=[],
    edges=[],
	escapes=[],
	coincidences=[],
    autoklms = false,
    xmin=0.0,
    xmax=missing,
    norm=:None,
    yscale=1.05,
    ytransform = identity,
	style=NeXLSpectrumStyle,
	palette=NeXLPalette
)::Plot

Required:

specs::AbstractVector{Spectrum};

Named:

klms = Union{Element, CharXRay}[ ]
edges = Union{Element, AtomicSubShell}[ ]
escapes = Union{Element, CharXRay}[ ],
coincidences = CharXRay[ ],
autoklms = false # Add KLMs based on elements in spectra
xmin = 0.0 # Min energy (eV)
xmax = missing # Max energy (eV) (defaults to max(:BeamEnergy))
norm = NoScaling() | ScaleDoseWidth() | ScaleDose() | ScaleSum() | ScaleROISum() | ScalePeak() | (<: SpectrumScaling)()
yscale = 1.05 # Fraction of max intensity for ymax over [max(lld,xmin):xmax]
ytransform = identity | log10 | sqrt | ??? # How to transform the counts data before plotting
style=NeXLSpectrumStyle (or another Gadfly.style)
palette = NeXLPalette | Colorant[ ... ] # Colors for spectra...
customlayers = Gadfly.Layer[] # Allows additional plot layers to be added

Plot a multiple spectra on a single plot using Gadfly.

source

These types define the different ways that spectra can be scaled when plotted using the Gadfly.plot(...) methods.

NeXLSpectrum.SpectrumScalingType

SpectrumScaling types are designed to rescale spectrum data primarily for plotting.

Implement

Base.show(io::IO, scn::SpectrumScaling)
scalefactor(sc::SpectrumScaling, spec::Spectrum)
source
NeXLSpectrum.ScaleWidthType

Scale to a constant dose⋅width (Counts/(nA⋅s/eV)) Requires a spectrum has both :ProbeCurrent & :LiveTime.

source
NeXLSpectrum.NeXLSpectrumStyleConstant

NeXLSpectrumStyle defines the default look-and-feel for Gadfly.plot(...) as applied to EDS spectra using the Gadfly.plot(...) functions implemented in NeXLSpectrum.

source

Spectrum Tabulation

NeXLUncertainties.asaMethod
NeXLUncertainties.asa(::Type{DataFrame}, spec::Spectrum; properties::Bool = false)

Converts the spectrum energy and counts data into a DataFrame.

source
NeXLUncertainties.asaMethod
NeXLUncertainties.asa(::Type{DataFrame}, spec::AbstractArray{Spectrum})::DataFrame

Returns a DataFrame that summarizes the list of spectra.

source

HyperSpectrum Manipulation

NeXLSpectrum.HyperSpectrumType
HyperSpectrum(arr::Array{T<:Real}, energy::EnergyScale, props::Array{Symbol, Any})

HyperSpectrum(energy::EnergyScale, props::Dict{Symbol,Any}, arr::Array{<:Real}; axisnames = ( "X", "Y", "Z", "A", "B", "C" ), fov = ( 1.0, 1.0, 1.0, 1.0, 1.0, 1.0 ), livetime=)
HyperSpectrum(energy::EnergyScale, props::Dict{Symbol,Any}, arr::AxisArray)
HyperSpectrum(energy::EnergyScale, props::Dict{Symbol,Any}, dims::NTuple{<:Integer}, depth::Int, type::Type{Real}; axisnames = ( "X", "Y", "Z", "A", "B", "C" ), fov = ( 1.0, 1.0, 1.0, 1.0, 1.0, 1.0 )

The HyperSpectrum struct represents a multi-dimensional array of Spectrum objects. The dimension of a HyperSpectrum may be 1 for a traverse or a line-scan, 2 for a spectrum image or >2 for time-series of spectrum images or multi-slice spectrum images.

The first constructor is used to create a HyperSpectrum from a raw Array of data. The second to construct a HyperSpectrum from another HyperSpectrum or an AxisArray. The third from a description of the intended contents.

  • axisnames: A list of the names by which the axis can be referred
  • fov: The full width of the dimension in mm.

HyperSpectra differ from Array{Spectrum} in that the spectra in a HyperSpectrum must share properties like :ProbeCurrent and :BeamEnergy. However, each pixel can have it's own livetime. HyperSpectrum objects can refer to line-scans (1D), spectrum images (2D), slice-and-view (3D), time sequenced images (3D), or higher dimension spectrum images.

Internally, HyperSpectrum reinterpretes an Array{T<:Real, N+1} as an Array{Spectrum{T<:Real},N-1}.

HyperSpectrum objects can be read from a RPL/RAW file (using readrplraw(filenamebase::AbstractString)) but can be constructed from any Array{<:Real}.

HyperSpectrum objects can be indexed using the standard Julia array idioms including a single integer index or a CartesianIndex. For example, to iterate over every spectrum in a HyperSpectrum

% Construct a 21 row × 19 column spectrum image with 2048 channels of [0,255].
hs = HyperSpectrum{es, props, (21,19), 2048, UInt8}
for idx in eachindex(hs)
    spec = hs[idx]   % get a Spectrum representing the 2048 channels of data at idx
    spec[22] = 1     % Set the 22nd channel to 1
end
% Indices into a HyperSpectrum are (row, column)
source
NeXLSpectrum.linescanFunction
linescan(hss::HyperSpectrum{T,2,3}, ci1::CartesianIndex{2}, ci2::CartesianIndex{2}, width::Int=1)

Extract pixels from hss along the line from ci1 to ci2 as a 1D HyperSpectrum. The width argument integrates the linescan along a line perpendicular to the primary axis. Only works on 2D SpectrumImages and for odd values of width. The algorithm does double count the occasional pixel but the length of each perpendicular is maintained at width. The AxisArrays.axes(...) scaling is maintained so lengths on the linescan can be compared to lengths on the original map.

source
NeXLSpectrum.blockFunction
block(hss::HyperSpectrum{T,N,NP}, steps::NTuple{N, Int})::HyperSpectrum{T,N,NP} where {T<:Real, N, NP}
block(hss::HyperSpectrum{T,N,NP}, step::Int) where {T<:Real, N, NP}

Reduce the size of a HyperSpectrum by summing together blocks of adjacent pixels. For example, steps=(4,4) would sum together blocks of 16 spectra in hss to form a single pixel in the resulting Hyperspectrum.

source
NeXLSpectrum.readrplrawFunction
readrplraw(rplfilename::AbstractString, rawfilename::AbstractString)
readrplraw(rplio::IO, rawio::IO)

Read a RPL/RAW file pair from IO or filename into an Array obect. The reader supports :bigendian, :littleendian ordering and :vector or :image alignment. Since the Array can be very large it is beneficial to release and collected the associated memory when you are done with the data by assigning nothing to the variable (or allowing it to go out of scope) and then calling GC.gc() to force a garbage collection.

* Ordering: The individual data items may be in `:littleendian` or `:bigendian`.
** `:littleendian` => Intel/AMD and others
** `:bigendian` => ARM/PowerPC/Motorola
* Alignment:  The data in the file can be organized as `:vector` or `:image`.  However, the data will be

reorganized into 'vector' format when returned as a Array. **:vector => Spectrum/data vector as contiguous blocks by pixel ** :image => Each channel of data organized in image planes * Data types: signed/unsigned 8/16/32-bit integers or 16-bit/32-bit/64-bit floats

readrplraw(filenamebase::AbstractString)::Array{<:Real}

Read the files filenamebase".rpl" and filenamebase".raw" into an Array. Maintains the data type of the values in the RAW file.

Standard LISPIX Parameters in .rpl File

.rpl files consist of 9 lines. Each line consists of a 'key'<tab>'value' where there is one and only one tab and possibly other space between the parameter name and parameter value. Parameter names are case-insensitive. The first line in the files is "key<tab>value". Subsequent lines contain the keys and values described in this table.

keyvaluedescription
width849pixels per row integer
height846rows integer
depth4096images or spec pts integer
offset0bytes to skip integer
data-length1bytes per pixel 1, 2, 4, or 8
data-typeunsignedsigned, unsigned, or float
byte-orderdont-carebig-endian, little-endian, or dont-care
record-byvectorimage, vector, or dont-care

This .rpl file indicates the image is 849 pixels wide and 846 pixels high, with 4096 levels in the depth dimension.

source
NeXLSpectrum.readptxFunction
readptx(
    fn::AbstractString, 
    scale::EnergyScale, # Energy scale for hyperspectrum
    props::Dict{Symbol,Any},
    nch::Int; # Number of channels in hyperspectrum
    blocksize = 1,  # Size of blocks to sum to create averaged hyperspectrum
    dets=[true,true,true,true], # Which detectors to include
    frames=1:typemax(Int)) # which frames to include in hyperspectrum
readptx(
    fn::AbstractString, 
    scale::EnergyScale, # Energy scale for hyperspectrum
    nch::Int; # Number of channels in hyperspectrum
    blocksize = 1,  # Size of blocks to sum to create averaged hyperspectrum
    dets=[true,true,true,true], # Which detectors to include
    frames=1:typemax(Int)) # which frames to include in hyperspectrum

Read a SEMantics .ptx file into a HyperSpectrum.

source
NeXLSpectrum.readhspyFunction
readhspy(filename::String, name=nothing)

Read a HyperSpectrum from a HyperSpy-style HDF5 file. The name argument allows the user to optionally specify an experiment to load.

source
NeXLSpectrum.planeFunction

plane(hss::HyperSpectrum, chs::AbstractUnitRange{<:Integer}, normalize=false)

Sums a contiguous range of data channels into an Array. The dimension of the result is one less than the dimension of the HyperSpectrum and is stored as a Float64 to ensure that not information is lost.

source

plane(hss::HyperSpectrum, ch::Int, normalize=false)

Extracts a single channel plane from a HyperSpectrum. The dimension of the result is one less than the dimension of the HyperSpectrum.

source
NeXLSpectrum.roiimageFunction
roiimage(hss::HyperSpectrum, chs::AbstractUnitRange{<:Integer})

Create a count map from the specified contiguous range of channels. (Accounts for differences in :LiveTime between pixels.)

source
roiimage(hss::HyperSpectrum, cxr::CharXRay)

Create a count map for the specified characteristic X-ray. By default, integrates for one FWHM at cxr. If hss[:Detector] is an EDSDetector, the FWHM is taken from it. Otherwise, a FWHM of 130 eV at Mn Kα is assumed.

source
NeXLSpectrum.compressFunction
compress(hss::HyperSpectrum)

Returns a HyperSpectrum with smaller or equal storage space to hss without losing or truncating any counts (note: AbstractFloat compresses to Float32 with loss of precision). Can change the storage type and/or reduce the depth of hss. If there is nothing that can be done to reduce the size of hss then compress(hss)===hss.

Example:

 hs = compress(hs) # Replace `hs` with a version that uses less memory (if possible).
source
NeXLSpectrum.maxpixelFunction
maxpixel(hss::HyperSpectrum, filt=ci->true)
maxpixel(hss::HyperSpectrum, cis::CartesianIndices, filt=ci->true)
maxpixel(hss::HyperSpectrum, mask::BitArray)
minpixel(hss::HyperSpectrum, filt=ci->true)
minpixel(hss::HyperSpectrum, cis::CartesianIndices, filt=ci->true)
minpixel(hss::HyperSpectrum, mask::BitArray)

Compute Bright's Max-Pixel derived signal for the entire HyperSpectrum or a rectanglar sub-region.

source
NeXLCore.colorizeFunction
colorize(hss::HyperSpectrum, cxrs::AbstractVector{CharXRay}, normalize=:All)

Create RGB colorized images from up to three CharXRay which define channels over which the count signal is integrated. normalize=:All puts the intensities on a common scale using the roiimages(...) method. Otherwise each image is scaled independently based on the brightest pixel.

source
NeXLSpectrum.labeledimagesFunction
labeledimages(labels::AbstractVector{<:AbstractString}, images::AbstractVector{<:AbstractArray}; ncols=3, halign = hleft)

Create a matrix of labeled images. Useful in Weave, Juno or Jupyter.

source
NeXLSpectrum.labeledimageFunction
labeledimage(label::String, image::Array)

Displays a string label below an image. Useful in Weave, Juno or Jupyter.

source
NeXLSpectrum.regionFunction
region(hss::HyperSpectrum{T, N}, ranges::AbstractRange...)::HyperSpectrum where {T<:Real,N}

Creates a view of a HyperSpectrum to represent the range of pixels within the hss HyperSpectrum. Does not copy the data or properties so any modifications to the region are also made to hss.

source
NeXLSpectrum.indexofmaxpixelFunction
indexofmaxpixel(hss::HyperSpectrum, ch::Int) # at channel `ch`
indexofmaxpixel(hss::HyperSpectrum) # all channels
indexofmaxpixel(hss::HyperSpectrum, ch::Int, cis::CartesianIndices)
indexofmaxpixel(hss::HyperSpectrum, cis::CartesianIndices)

Find the indices producing the maximum value in data[ch] or data[:] within 'cis' or full spatial dimensions.

source
NeXLSpectrum.roiimagesFunction
roiimages(hss::HyperSpectrum, achs::AbstractVector{<:AbstractUnitRange{<:Integer}})

Create an array of Gray images representing the intensity in each range of channels in in achs. They are normalized such the the most intense pixel in any of them defines white. (Accounts for differences in :LiveTime between pixels.)

source
roiimages(hss::HyperSpectrum, cxrs::AbstractVector{CharXRay})

Create an array of Gray images representing the intensity in each of the CharXRay lines in cxrs. They are normalized such the the most intense pixel in any of them defines white. By default, integrates for one FWHM at cxr. If hss[:Detector] is an EDSDetector, the FWHM is taken from it. Otherwise, a FWHM of 130 eV at Mn Kα is assumed.

source
NeXLSpectrum.livetime!Function
livetime!(hss::HyperSpectrum, lt::AbstractFloat, idx...)
livetime!(hss::HyperSpectrum{T,N}, lt::AbstractFloat) # All pixels to lt

Set the livetime on a per pixel basis.

source

Fitting Filter

Core methods for constructing FilterFitPackets and fitting spectra.

NeXLSpectrum.referenceFunction
reference( elm::Element, spec::Spectrum, mat::Material=spec[:Composition]; pc = nothing, lt = nothing, e0 = nothing, coating = nothing)::ReferencePacket
reference(els::AbstractVector{Element}, spec::Spectrum, mat::Material = spec[:Composition]; pc = nothing, lt = nothing, e0 = nothing, coating = nothing)::Vector{ReferencePacket}

Construct a ReferencePacket from a Spectrum collected from the specified Material for the specified Element. Often used with references(...) to build FilterFitPackets.

Optional named arguments pc, lt, e0, coating allow you to specify the probe current, live time, beam energy and sample coating.

source
NeXLSpectrum.referencesFunction
references(refs::AbstractVector{ReferencePacket}, det::EDSDetector)::FilterFitPacket
references(refs::AbstractVector{ReferencePacket}, fwhm::Float64)::FilterFitPacket

Constructs a FilterFitPacket from a vector of ReferencePackets. Each ReferencePacket represents a single ROI for an element. It is possible more than one ReferencePacket might be defined for an elemental ROI. In this case, the ReferencePacket with the lower index will take preference over later ones. This allows you to fill in only the missing elemental ROIs using spectra collected from alternative materials. For example, a spectrum from F₂Fe is suitable for the Fe K-lines but not the Fe L-lines. So we might specify F₂Fe first to specify the references for the Fe K-lines and then fill in the L-lines with a spectrum from pure Fe.

source
references(refs::Vector{DirectRefInit}, det::Detector, resp::Matrix{Float64}; minE=0.5e3 )

Use along with direct(...) to build a collection of direct fitting references.

Example:

> detu = matching(unk, 132.0, 110)
> resp = detectorresponse(detu, SDDEfficiency(ModeledWindow(MoxtekAP33())))
> drefs = references([ 
      direct(n"O", stds[1], mat"Al2O3"),
      direct(n"Al", stds[1], mat"Al2O3"),
      direct(n"Ba", stds[2], mat"BaF2"),
      direct(n"Ca", stds[3], mat"CaF2"),
      direct(n"Fe", stds[4], mat"Fe"),
      direct(n"Si", stds[5], mat"Si") ],
      detu, resp
  )
source
NeXLSpectrum.FilterFitPacketType

Represents the processed spectral data necessary to efficiently filter-fit one or more unknown spectra. A FilterFitPacket{S<:Detector, T<:AbstractFloat} contains the data necessary to filter the unknown and to apply pre-filtered references. The type T <: AbstractFloat allows you to specify the bit-depth used to perform subsequent calculations using this FilterFitPacket. Float32 uses less memory and is a little faster than Float64. If there are duplicate FilteredReference for an elemental ROI, the preference is for the first one. This allows you to fill in unavailable "FilteredReference" elemental ROIs with more general ones.

source
NeXLSpectrum.suitabilityFunction
suitability(elm::Element, mats::AbstractSet{<:Material}, det::Detector; maxE=30.0e3)
suitability(elm::Element, det::Detector; maxE=30.0e3, minC=0.1)

Tabulates the characteristic X-ray peaks for the Element for which there are suitable materials in mats for the specified detector. The second form uses a default set of Materials in the file NeXLCore "standards.txt".

This function is helpful for determining which Materials are suitable to act as fitting standards for the specified Element. It shows how NeXLSpectrum will break up the characteristic peaks associated with elm into contiguous regions each of which will be fit independently. NeXLSpectrum attempts to break each element into as many independent regions as possible dependent on the resolution of the specified EDSDetector. If there is an interference between one of the other elements in the Material and elm then this peak will not be suitable as a fitting standard. However, it can be used as a similar standard.

source
NeXLSpectrum.suitableforFunction
suitablefor(
    elm::Element, 
    matOrElms::Union{Material, AbstractSet{Element}}, 
    det::Detector; 
    maxE::Float64 = 1.0e6, 
    ampl::Float64 = 1.0e-5,
    warnme::Bool = true
)::Vector{Tuple{Vector{CharXRay}, UnitRange{Int64}}}

Given a material or collection of Element which ROIs for element elm is the material a suitable reference on the detector det? This function provides informational, warning and error messages depending upon the suitability of the material.

source
NeXLSpectrum.fit_spectrumFunction
fit_spectrum(
    hs::Spectrum|HyperSpectrum,
    vq::VectorQuant{S <: Detector, T <: AbstractFloat},
    zero = x -> max(Base.zero(T), x)
)

Fit the spectrum or hyper-spectrum using the vector-quant algorithm. The function zero is applied to the resultant k-ratios before they are returned. The default simply sets negative k-ratios to 0.0. zero=identity would leave the negative k-ratios as such.

source
fit_spectrum(unk::Spectrum{T}, drefs::DirectReferences)::DirectFitResult{T} where { T <: Real }

Fits the direct references to the unknown spectrum returning a DirectFitResult struct containing k-ratios and other output quantities from the fit.

Surprisingly, this function requires that unk[:Composition] is defined as a Material with an estimate of the composition of the material from which unk was collected. This information is necessary to model the continuum background. Yes, this seems circular (and it is.) One option is to perform a filter-fit first, quantify the filter fit and then use this as your estimate.

source
NeXLSpectrum.FilteredReferenceType
FilteredReference

Represents the filtered reference spectrum over an ROI. Carries the minimal data necessary to support filter-fitting a single region-of-interest (continguous range of channles) and computing useful output statistics.

source
NeXLCore.elmsFunction
elms(spec::Spectrum, withcoating = false)::Set{Element}

Returns a set of the elements associated with this spectrum. withcoating determines whether the coating elements are also added.

source

A list of elements for which there are filtered references.

source
NeXLSpectrum.missingReferencesFunction
missingReferences(ffp::FilterFitPacket, elms::Vector{Element}, e0::Float64, ampl=1.0e-5)

Returns a Vector{Tuple{Vector{CharXRay}, UnitRange{Int64}}} containing the ROIs for which a FilteredReference is missing from the FilterFitPacket.

source
NeXLSpectrum.FilterFitResultType
FilterFitResult

Represents the result of fitting a FilteredUnknownW to a FilteredUnknown.

Struct elements

label::UnknownLabel  # Identifies the unknown
kratios::UncertainValues # Labeled with ReferenceLabel objects
roi::UnitRange{Int} # Range of channels fit
raw::Vector{T} # Raw spectrum data
residual::Deferred # Residual spectrum
peakback::Deferred # {Dict{ReferenceLabel,NTuple{3,Float64}}} with peak counts, background counts and counts/(nAs)

Use asa(DataFrame, ffr::FilterFitResult) to summarize in tabular form.

source
NeXLSpectrum.kratiosFunction
kratios(ffr::FitResult)::Vector{KRatio}

The k-ratios associated with each CharXRayLabel as a vector 'KRatio' objects.

source
NeXLSpectrum.residualFunction
residual(ffr::FilterFitResult)::Spectrum

A Spectrum containing the histogram representing the unknown spectrum minus the fitted characteristic peaks shapes times the best fit coefficient.

source
NeXLSpectrum.filteredresidualFunction
filteredresidual(fit::FilterFitResult, unk::FilteredUnknown, ffs::AbstractVector{FilteredReference})::Vector{Float64}

Computes the difference between the best fit and the unknown filtered spectral data.

source
NeXLUncertainties.extractFunction
extract(fd::FilteredReference{T}, roi::UnitRange{Int})::Vector{T} where { T <: AbstractFloat }
extract(fd::FilteredUnknown, roi::UnitRange{Int})::Vector{T} where { T <: AbstractFloat }

Extract the filtered data representing the specified range. roi must fully encompass the filtered data in fd.

source
NeXLSpectrum.fit_spectraFunction
fit_spectrum(spec::Spectrum, ffp::FilterFitPacket)::FilterFitResult
fit_spectrum(specs::AbstractVector{<:Spectrum}, ffp::FilterFitPacket)::Vector{FilterFitResult}
fit_spectra(specs::AbstractVector{<:Spectrum}, ffp::FilterFitPacket)::Vector{FilterFitResult}

Fit a Spectrum or a vector of Spectrum using the specified FilterFitPacket. The result is a FilterFitResult structure which contains k-ratios, residuals, etc.

fit_spectrum(hs::HyperSpectrum, ffp::FilterFitPacket; mode::Symbol=:Fast, zero = x -> max(0.0, x))::Array{KRatios}
  • mode = :Fast - Uses precomputed, filtered "vector" fit method. No uncertainties are available.
  • mode = :Intermediate - Uses an optimized full fit without refits for negative k-ratios.
  • mode = :Full - Uses the full single spectrum fit algorithm including refitting when one or more k-ratio is less than zero.

Performs a filtered fit on a hyperspectrum returning an Array{KRatios}.

Selecting a mode: :Fast is good for generating k-ratio maps or exploratory analysis of a k-ratio map. :Full is best when a quantitative map of a high count hyperspectrum is desired. Fit results for major elements are similar for all three but differ for minor and trace elements. Particularly when a k-ratio is slightly negative. This negative k-ratio can effect the other k-ratios. :Fast also works less well when many elements (>>10) (particularly interfering elements) are included in the fit. Unfortunately, :Intermediate and :Full slow down when many elements are fit - O(n(elements)²).

The following timing on a 512 x 512 x 2048 hyperspectrum fitting 21 elements with 33 distinct ROIs on a fast 6-core laptop with 16 GiB memory give a relative feel for the speed of each algorithm. Yes, :Fast is approximately 30 × faster than :Intermediate and used almost 80x less memory. (64-bit references, 32-bit references use about 1/2 the memory and take about 4/5 the time.)

mode=ThreadsRun time (s)AllocationsAllocated (GiB)GC time
:Fast611.42.11 M4.723.0%
:Intermediate6342.713.11 M364.44.2%
:Full6574.66.2 G862.95.5%
:Fast125.52.1 M4.721.8%
:Intermediate11064.613.1 M364.40.9%
:Full12186.26.2 G862.12.8%
source

Filter Fit Tabulation

NeXLUncertainties.asaMethod
NeXLUncertainties.asa(::Type{DataFrame}, ffp::FilterFitPacket)

Summarize the FilteredReference structs within a FilterFitPacket as a DataFrame.

source
NeXLUncertainties.asaMethod
NeXLUncertainties.asa(::Type{DataFrame}, ffrs::AbstractVector{<:FitResult}; charOnly = true, withUnc = false, format = :normal # :pivot or :long)

Return generic FitResult as a DataFrame.

Format:

  • :normal - One row per spectrum, one column per k-ratio
  • :pivot - One row per ROI, one column per spectrum (optional: column for 1σ uncertainty on k-ratio)
  • :long - One row per spectrum, feature, measured k-ratio
source
NeXLUncertainties.asaMethod
NeXLUncertainties.asa(
    ::Type{DataFrame}, 
    ffr::FilterFitResult; 
    charOnly::Bool=true, 
    material=nothing,
    mc = XPP, fc=ReedFluorescence,
    columns = ( :counts, ) # Selected from ( :roi, :peakback, :counts, :dose )
)::DataFrame

Tabulate details about each region-of-interest in the 'FilterFitResult' in a 'DataFrame'.

  • If charOnly then only display characteristic X-ray data (not escapes etc.)
  • If material is a Material then the computed k-ratio (KCalc) will also be tabulated along with kmeas/kcalc (KoKCalc).
  • columns - Select optional column outputs (see below)

Columns:

  • Spectrum : UnknownLabel - Identifies the fit spectrum
  • Feature : Label - Identifies the fit feature (Vector{CharXRay} or other)
  • Reference: String - Name of reference against which :Spectrum was fit over :Feature
  • K : Float64 - The multiplicative fit constant
  • dK : Float64 - The 1σ uncertainty in :K
  • :Start : Int - Start index for fit channels (:roi ∈ columns)
  • :Stop : Int - Stop index for fit channels (:roi ∈ columns)
  • :Counts : Float64 - Total counts in characteristic peak (peak-back) (:peakback || :counts ∈ columns)
  • :Back : Float64 - Total counts in background under the characteristic peak (:peakback ∈ columns)
  • :PtoB : Float64 - Peak-to-Background assuming 10 eV/channel (:peakback ∈ columns)
  • :KCalc : Float64 - Computed k-ratio assuming a composition. (Requires material argument to be specified.)
  • :KoKcalc : Float64 - Ratio of measured/computed k-ratio. (Requires material argument to be specified.)
  • :LiveTime : Float64 - Acquisiton live time (s) (:dose ∈ columns)
  • :ProbeCurrent : Float64 - Probe current (nA) (:dose ∈ columns)
  • :DeadPct : Float64 - Dead time in ProbeCurrent (:dose ∈ columns)
  • :RefCountsPernAs : Float64 - Estimated counts in :Reference in :Feature per unit dose. (:counts ∈ columns)
  • :CountsPernAs : Float64 - Estimated counts in :Spectrum in :Feature per unit dose. (:counts ∈ columns)
source
NeXLUncertainties.asaMethod
NeXLUncertainties.asa(::Type{DataFrame}, fr::FitResult; withUnc = false)

Summarize the ROI and k-ratio data within a FitResult structure as a DataFrame.

source

Filter Fit Plotting

Gadfly.plotMethod
Gadfly.plot(vq::VectorQuant, chs::UnitRange)

Plots the "vectors" used to quantify various elements/regions-of-interest over the range of channels specified.

source
Gadfly.plotMethod
plot(ff::TopHatFilter, fr::FilteredReference)

Plot a color map showing the filter data relevant to filtering the specified FilteredReference.

source
Gadfly.plotMethod
Gadfly.plot(fr::FilteredReference; palette = NeXLPalette))

Plot a filtered reference spectrum.

source
Gadfly.plotMethod
Gadfly.plot(ffp::FilterFitPacket; kwargs...)

Plots the reference spectra which were used to construct a FilterFitPacket.

source
Gadfly.plotMethod
Gadfly.plot(
    ffr::FilterFitResult,
    roi::Union{Nothing,AbstractUnitRange{<:Integer}} = nothing;
    palette = NeXLPalette,
    style = NeXLSpectrumStyle,
    xmax::Union{AbstractFloat, Nothing} = nothing,
    comp::Union{Material, Nothing} = nothing,
    det::Union{EDSDetector, Nothing} = nothing,
    resp::Union{AbstractArray{<:AbstractFloat,2},Nothing} = nothing,
    yscale = 1.0
)

Plot the sample spectrum, the residual and fit regions-of-interests and the associated k-ratios.

source

Advanced Filter Fitting

NeXLSpectrum.TopHatFilterType

The TopHatFilter{T <: AbstractFloat} struct represents a zero-sum symmetric second-derivative-like filter that when applied to spectral data has the property of suppressing constant and slowly varying signals (like the continuum) while retaining a linear signal for faster changing signals like the characteristic peaks.

See

  • F. H. Schamber Proc Symposium of "X-ray Fluorscence Analysis on Environmental Samples" Chapel Hill 1976 T Dzubay Ed.
  • P. Statham Anal Chem 49 no 14 Dec 1977

The TopHatFilter struct optimizes the memory and CPU use when applying top-hat filters to spectrum data.

The easiest way to implement a top-hat filter is as matrix F. The rows represent the filters. The product of the filter and the data vector is the filtered spectrum. The product of the filter times a diagnonal matrix constructed from the data times the transpose of the filter is the covariance of the filtered data. The diagonal matrix constructed from the spectrum data is the covariance matrix associated with the spectrum data because the channels in the spectrum data are independent (thus the matrix is diagnonal) and the magnitude equals the counts in each channels because the spectrum data is nominally Poissonian and in the large number limit, the variance of a Poissonian random variable is the number itself (σ=sqrt(N) => Var = N)

Notes on memory and code optimization: The filter matrix is banded diagonal. Approximately, 2.5% of the elements are non-zero. This suggest use of the BandedMatrix type. The most expensive operation is calculating F⋅D⋅Fᵀ, the covariance matrix of the filtered data. D is a diagonal matrix and so computing each element in F⋅D⋅Fᵀ reduces to a sum over a single variable. Furthermore, the weighted least squares fit doesn't require the full F⋅D⋅Fᵀ, just diag(F⋅D⋅Fᵀ). However, it turns out that we can do better implementing our own banded matrix type largely because D is fully diagonal and the matrix product F⋅D⋅Fᵀ reduces down to a sum over a single variable. The product F⋅d and F⋅D⋅Fᵀ are readily implemented as element-by-element multiplies and sums. Thus storing the filter as offsets and row filters is efficient in both memory and CPU use.

source
NeXLSpectrum.GaussianFilterType
GaussianFilter

A Gaussian-shaped filter that varies in width with the FWHM of the detector. The Gaussian is offset to ensure the sum of the filter elements is zero.

source
NeXLSpectrum.tophatfilterFunction
tophatfilter(
    charLabel::Union{CharXRayLabel,EscapeLabel, ReferenceLabel}
    filt::TopHatFilter{T},
    scale::T = one(T),
    tol::T = (T==Float32 ? 1.0f-5 : 1.0e-6)
)::FilteredReference

For filtering an ROI on a reference spectrum. Process a portion of a Spectrum with the specified filter. Use a simple edge-based background model.

tophatfilter(
    labels::AbstractVector{ReferenceLabel},
    filt::TopHatFilter{T},
    scale::T = one(T),
    tol::T = (T==Float32 ? 1.0f-5 : 1.0e-6)
)::Vector{FilteredReference{T}}
source
tophatfilter(spec::Spectrum, filt::TopHatFilter{T}, scale::T = one(T))::FilteredUnknown where { T <: AbstractFloat }

For filtering the unknown spectrum. Defaults to the weighted fitting model.

source
tophatfilter(::Type{FilteredUnknownW}, spec::Spectrum, thf::TopHatFilter, scale::Float64=1.0, tol::Float64 = 1.0e-4)::FilteredUnknownW

For filtering the unknown spectrum. Process the full Spectrum with the specified filter for use with the weighted least squares model.

source
NeXLSpectrum.buildfilterFunction
buildfilter(
    [ ::Type{T} = Float64],
    ty::Type{<:TopHatFilterType},
    det::Detector,
    a::AbstractFloat = 1.0, # Top
    b::AbstractFloat = 2.0,  # Base
)
buildfilter(::Type{T}, det::Detector, a::AbstractFloat = 1.0, b::AbstractFloat = 2.0) where { T<:AbstractFloat }
buildfilter(det::Detector, a::AbstractFloat = 1.0, b::AbstractFloat = 2.0)::TopHatFilter{Float64}
buildfilter(ty::Type{<:TopHatFilterType}, det::Detector, a::AbstractFloat = 1.0, b::AbstractFloat = 2.0)::TopHatFilter{Float64}

Build a top-hat filter for the specified detector with the specified top and base parameters. The default element type is Float64 and the default shape model (TopHatFilterType) is VariableWidthFilter.

source
buildfilter(::Type{GaussianFilter}, det::Detector, a::AbstractFloat=1.0, b::AbstractFloat=6.0)::TopHatFilter

Build a top-hat filter with Gaussian shape whose width varies with the detector's resolution as a function of X-ray energy for the specified detector with the specified top and base parameters. The a parameter corresponds to the filter width relative to the detector resolution expressed as Gaussian width. So a=1 is a filter whose width equals the detector resolution at each energy. The b parameter is the extent of the filter in Gaussian widths. The default a=1, b=5 corresponds to a filter that has the same resolution as the detector and an extent of 2.5 Gaussian widths above and below the center channel.

source
NeXLSpectrum.FilteredUnknownWType
FilteredUnknownW

Represents the unknown in a filter fit using the weighted fitting model. This is an approximation that produces over optimistic resulting covariance matrix.

source
NeXLSpectrum.filterfitFunction
filterfit(unk::FilteredUnknownW, ffs::AbstractVector{FilteredReference}, forcezeros = true)::FilterFitResult

Filter fit the unknown against ffs, an array of FilteredReference and return the result as an FilterFitResult object. By default use the generalized LLSQ fitting (pseudo-inverse implementation).

This function is designed to reperform the fit if one or more k-ratio is less-than-or-equal-to zero. The FilteredReference corresponding to the negative value is removed from the fit and the fit is reperformed. How the non-positive value is handled is determine by forcezeros. If forcezeros=true, then the returned k-ratio for the non-positive value will be set to zero (but the uncertainty remains the fitted one). However, if forcezeros=false, then the final non-positive k-ratio is returned along with the associated uncertainty. forcezeros=false is better when a number of fit k-ratio sets are combined to produce an averaged k-ratio with reduced uncertainty. forcezeros=true would bias the result positive.

source
NeXLSpectrum.isvisibleFunction
isvisible(cxrs::AbstractVector{<:SpectrumFeature}, det::Detector)

Returns the characteristic x-rays that are visible on the specified detector (ie. Between the LLD and the maximum channel).

source
isvisible(sf::SpectrumFeature, det::Detector)

Is sf visible on the specified Detector?

source
NeXLSpectrum.ReferenceLabelType
ReferenceLabel

A label associated with reference spectra. The label encapsulates the original spectrum and the range of channels represented by this reference object. structs that extend ReferenceLabel should have .roi, .spectrum and ".hash" fields.

source
NeXLSpectrum.SpectrumFeatureType
SpectrumFeature

A union representing the different type of peak-like features (helpful and harmful) that can appear in a spectrum.

source
NeXLSpectrum.CharXRayLabelType
CharXRayLabel

A ReferenceLabel that represents a reference spectrum or reference properties associated with a set of characteristic x-rays (CharXRay) objects over a contiguous range of spectrum channels.

source
NeXLSpectrum.EscapeLabelType
EscapeLabel

A ReferenceLabel<:FilteredLabel that Represents a reference spectrum associated with an escape peak from a set of characteristic x-rays (CharXRay) objects over a contiguous range of spectrum channels.

source
NeXLSpectrum.charXRayLabelsFunction
charXRayLabels(#
  spec::Spectrum, #
  elm::Element, #
  allElms::AbstractSet{Element}, #
  det::Detector, #
  ampl::Float64, #
  maxE::Float64=1.0e6)::Vector{SpectrumFeature}

Creates a vector of CharXRayLabel objects associated with 'elm' for a spectrum containing the elements 'allElms' assuming that it was collected on 'det'. ROIs in which other elements from 'allElms' interfere with 'elm' will not be included.

source
NeXLSpectrum.directFunction
direct(elm::Element, spec::Spectrum, mat::Material=spec[:Composition])
direct(elm::Element, specfile::String, mat::Material)
direct(elm::Element, specfile::String)

Construct a struct to represent a direct-fit reference. Use with references(...) to construct a reference set which may be used to fit multiple references to an unknown spectrum using a "direct" linear fit.

source
NeXLSpectrum.detectFunction
detect(emitted::Dict{<:XRay,<:Real}, det::EDSDetector, response::Matrix{Float64}; noise=false)

Returns a Spectrum as though the intensities in emitted were detected on the specified detector. noise=true will add Poisson noise to the resulting measured Spectrum.

source
NeXLSpectrum.filterreferenceFunction
filterreference(
    filt::TopHatFilter{T},
    spec::Spectrum,
    elm::Element,
    allElms::[AbstractSet{Element}|Material]
    props::Dict{Symbol,<:Any} = Dict{Symbol,Any}(),
    withEsc::Bool = false,
)
filterreferences(
    filt::TopHatFilter,
    refs::Tuple{Spectrum,Element,Material}...;
    props::Dict{Symbol,<:Any} = Dict{Symbol,Any}(),
    withEsc::Bool = false,
)
source

Matrix Correction

NeXLMatrixCorrection.quantifyFunction
NeXLMatrixCorrection.quantify(
    ffr::FitResult,
    iteration::Iteration = Iteration(mc=XPP, fc=ReedFluorescence, cc=Coating);
    strip::AbstractVector{Element} = Element[],
    kro::KRatioOptimizer = SimpleKRatioOptimizer(1.5),
    coating::Union{Nothing, Pair{CharXRay, <:Material}}=nothing
)::IterationResult

Facilitates converting FilterFitResult or BasicFitResult objects into estimates of composition by extracting k-ratios from measured spectra and applying matrix correction algorithms.

source
NeXLMatrixCorrection.quantify(
    spec::Union{Spectrum,AbstractVector{<:Spectrum}},
    ffp::FilterFitPacket;
    strip::AbstractVector{Element} = Element[],
    iteration::Iteration = Iteration(mc=XPP, fc=ReedFluorescence, cc=Coating),
    kro::KRatioOptimizer = SimpleKRatioOptimizer(1.5),
    coating::Union{Nothing, Pair{CharXRay, <:Material}}=nothing
)::IterationResult

Failitates quantifying spectra. First filter fits and then matrix corrects.

source
NeXLMatrixCorrection.estimatecoatingFunction
NeXLMatrixCorrection.estimatecoating(fr::FitResult, substrate::Material, coating::Material, cxr::CharXRay, mc::Type{<:MatrixCorrection}=XPP)::Film

Estimate the mass-thickness of coating on bulk substrate using the X-ray cxr for which there is KRatio in fr. The result is assigned to the properties of fr so that it can be account for in the matrix correction process.

This method is intended for use on standards where the substrate composition is known a priori.

source

Standardization

NeXLSpectrum.extractStandardsFunction
extractStandards(ffr::FitResult, elm::Element, mat::Material)::Vector{KRatio}

Extract a Vector{KRatio} for elm::Element from a ffr::FilterFitResult measured from mat::Material.

source
extractStandards(ffr::FitResult, cxrs::AbstractVector{CharXRay}, mat::Material)::Vector{KRatio}

Extract a KRatio for the CharXRay from a FilterFitResult associated with the Material.

source
NeXLCore.standardizeFunction
NeXLCore.standardize(ffr::FilterFitResult{T}, standard::FilterFitResult{T}, material::Material, els=elms(material))::FilterFitResult{T}
NeXLCore.standardize(ffrs::Vector{FilterFitResult{T}}, standard::FilterFitResult{T}, material::Material, els=elms(material))::Vector{FilterFitResult{T}}
NeXLCore.standardize(ffr::FilterFitResult{T}, standards::AbstractArray{KRatio})::FilterFitResult{T}
NeXLCore.standardize(ffrs::Vector{FilterFitResult{T}}, standards::AbstractArray{KRatio})::Vector{FilterFitResult{T}}

Apply the standard KRatios to the FilterFitResult producing a re-standardized FilterFitResult.

source

EDS Detectors

NeXLSpectrum.DetectorType
Detector

An abstract type defining the characteristics of an X-ray detector.

Implements:

channelcount(det::Detector)::Int
scale(det::Detector)::EnergyScale
resolution(eV::Float64, det::Detector)::Float64 # FWHM at eV
energy(ch::Int, det::Detector)::Float64
channel(eV::Float64, det::Detector)::Int
profile(energy::Float64, xrayE::Float64, det::Detector)
lld(det::Detector)::Int
isvisible(sf::SpectrumFeature, det::Detector)
source
NeXLSpectrum.EDSDetectorType
EDSDetector

Types extending EDSDetector must have member variables

channelcount::Int # Number of channels
scale::EnergyScale # Detector calibration funtion
resolution::Resolution # Detector lineshape function
lld::Int # low level discriminator
source
NeXLSpectrum.resolutionFunction
resolution(eV::Float64, res::Resolution)
resolution(eV::Float64, det::EDSDetector)

The FWHM at eV for the <:Resolution model.

source
NeXLSpectrum.ResolutionType
Resolution

An abstract type describing the channel dependence of the resolution of an EDS detector.

Implements:

resolution(eV::Float64, res::Resolution)::Float # Resolution at specified energy
profile(energy::Float64, xrayE::Float64, res::Resolution) # Amplitude for a signal at the specified energy at the specified energy
extent(xrayE::Float64, res::Resolution, ampl::Float64)::Tuple{2,Float} # The range of channels over which the signal exceeds ampl
source
NeXLSpectrum.simpleEDSwICCFunction
simpleEDSwICC(chCount::Integer, width::Float64, offset::Float64, fwhmatmnka::Float64, lld::Int=channel(150.0 eV))

Construct simple model of an EDS detector with incomplete charge collection at low X-ray energies.

source
NeXLSpectrum.matchesFunction
matches(spec::Spectrum, det::Detector, tol::Float64 = 1.0)::Bool
matches(spec1::Spectrum, spec2::Spectrum, tol::Float64 = 1.0)::Bool

Does the calibration of the Spectrum (approximately) match the calibration of the Detector or the other Spectrum?

source
matches(cxrl::CharXRayLabel, std::KRatio)
matches(cxrl::CharXRayLabel, std::StandardLabel)

Does the measured CharXRayLabel match the standard?

source
NeXLSpectrum.TabulatedWindowType
TabulatedWindow(wt::WindowType)

Construct a model of the window transmission based on vendor-supplied tabulations of transparency. At the end of the user supplied data, the transmission function is extended by matching the ModeledWindow with the tabulation.

source
NeXLSpectrum.DirectReferenceType

A DirectReference represents the data extracted from a reference spectrum necessary to perform a "direct" (background-eliminated direct fit).

source
NeXLSpectrum.AmptekC2Type

Create modeled windows for the Ametek C2 Si₃N₄ windows according to the specifications here: https://www.amptek.com/products/accessories-for-xrf-eds/c-series-low-energy-x-ray-windows#Specifications

Light transparent window often used for SEM detectors.

source
NeXLSpectrum.WindowTypeType

WindowType is the abstract type for models or types of X-ray windows. These types distinguish between the window type and the implementation. Often, there are both tabulated transmission functions from the vendor and calculated transmission functions based on the construction of the window. Use the TabulatedWindow or ModeledWindow to instantiate AbstractWindow which implements the transmission(wnd::AbstractWindow, energy::Float64, angle::Float64 = π / 2) method.

Predefined WindowTypes are MoxtekAP33, MoxtekAP5, AmptekC1, AmptekC2, BerylliumWindow.

source
NeXLSpectrum.DirectReferencesType

A DirectReferences is a packet of DirectReference and detector information. It can be used to fit the corresponding peaks in an "unknown" spectrum.

source
NeXLSpectrum.AbstractWindowType

AbstractWindow is the base type for TabulatedWindow, ModeledWindow and NoWindow. You can view the window transmission function using Gadfly and the method:

plot(wind::Union{AbstractWindow, AbstractArray{<:AbstractWindow}}; xmax=20.0e3, angle=π/2, style=NeXLSpectrumStyle)

Window types are identified by types based on WindowType like MoxtekAP33, MoxtekAP5, AmptekC1, AmptekC2, BerylliumWindow.

An implementation of an AbstractWindow is instantiate using code like TabulatedWindow(MoxtekAP33()) or ModeledWindow(AmptekC1()). In addition, there is a NoWindow() type that implements a 100% transparent window.

AbstractWindow types primarily implement NeXLCore.transmission(wnd::AbstractWindow, energy::Float64, angle::Float64 = π / 2)

source
NeXLSpectrum.ModeledWindowType
ModeledWindow(wt::WindowType)

This type models a window using the materials and thicknesses provided by the vendor. If accomodates grids by a simplistic mechanism of assuming an open area.

source
NeXLSpectrum.AmptekC1Type

Create modeled windows for the Ametek C1 Si₃N₄ windows according to the specifications here: https://www.amptek.com/products/accessories-for-xrf-eds/c-series-low-energy-x-ray-windows#Specifications

Light-tight (solar-blind) window. Often used for environmental XRF units.

source

EDS Detector Plotting

Gadfly.plotMethod
Gadfly.plot(deteff::Union{DetectorEfficiency,AbstractVector{DetectorEfficiency}}, emax = 20.0e3)

Plots the detector efficiency function assuming the detector is perpendicular to the incident X-rays.

source
Gadfly.plotMethod
Gadfly.plot(deteff::Union{DetectorEfficiency,AbstractVector{DetectorEfficiency}}, emax = 20.0e3)

Plots the detector efficiency function assuming the detector is perpendicular to the incident X-rays.

source

Energy Axis Scales for EDS Detectors

Structures and functions that implement the energy scale functions for EDS detectors.

NeXLSpectrum.EnergyScaleType
EnergyScale

An EnergyScale is a way of representing the energy axis associated with X-ray data. The scale may be linear, polynomial or ??? to handle the various different non-linearities that happen with EDS detectors plus we can also handle WDS wavescans.

Implements:

channel(::Type{Float64}, eV::AbstractFloat, sc::EnergyScale)::Float64
energy(ch::Int, sc::EnergyScale)::Float64
source
NeXLSpectrum.fwhmFunction
fwhm(gauss::Float64)

Converts Gaussian width to full-width half-max. See also gaussianwidth

source

Multi-Detector Functions

Functions for interpreting and manipulating spectr collected on multiple detectors simultaneously.

NeXLSpectrum.multicompareFunction
multicompare(specs::AbstractArray{Spectrum{T}}) where {T <: Real}

Compares the intensity for the spectra in specs against the mean intensity on a channel-by-channel basis. Compute the ratio for each channel in each spectrum of the spectrum intensity over the mean intensity for that channel. You expect the ratio to be unity when the spectra are identical and deviate from unity when the spectra are different.

source
NeXLSpectrum.multimeanFunction
multimean(specs::Spectrum{T})::Spectrum{T} where {T <: Real}

Average on a channel-by-channel basis spectra collected from multiple detectors simultaneously from the same electrons interacting with the same material for the real-time. The detectors should be calibrated close to identically to maintain the detector resolution and peak positions.

source
NeXLSpectrum.multisumFunction
multisum(specs::Spectrum{T})::Spectrum{T} where {T <: Real}

Sum together spectra collected from multiple detectors simultaneously from the same electrons interacting with the same material for the real-time. The detectors should be calibrated close to identically to maintain the detector resolution and peak positions.

source
NeXLSpectrum.plot_multicompareFunction
plot_multicompare(specs::AbstractArray{Spectrum{T}}; minE=200.0, maxE=0.5*specs[1][:BeamEnergy]) where { T<: Real}

Compare spectra collected simultaneously on multiple detectors in a single acquisition.

source

Bremsstrahlung

NeXLSpectrum.continuumroisFunction
continuumrois(elems, det::EDSDetector, minE::Float64, maxE::Float64)

Compute the ROIs for the contiguous continuum regions for the specified elements elems on an EDSDetector for the specified range of energies.

source
NeXLSpectrum.generatedFunction
generated(cm::ContinuumModel, ea::Float64)

Compute the intensity of the measured continuum generated from the material and conditions specified in the continuum model object at the specified measured energy ea.

source
NeXLSpectrum.emittedFunction
emitted(cm::ContinuumModel, ea::Float64)

Compute the intensity of the measured continuum emitted from the material and conditions specified in the continuum model object at the specified measured energy ea.

source
NeXLSpectrum.fitcontinuumFunction
fitcontinuum(
  spec::Spectrum,
  resp::AbstractArray,
  rois::Vector{UnitRange};
  brem::Type{<:NeXLBremsstrahlung} = Castellano2004a,
  mc::Type{<:MatricCorrection} = Riveros1993,
)

Fit a continuum model to the specified range of channels (`rois`).  The `resp` argument is a matrix which describes

the detector response on a channel-by-channel basis. It can be calculated from an EDSDetector and an DetectorEfficiency using resp = NeXLSpectrum.detectorresponse(det, eff). The Spectrum object must have the :Composition, :BeamEnergy and :TakeOffAngle properties defined.

source
fitcontinuum(
  spec::Spectrum,
  det::EDSDetector,
  resp::AbstractArray{<:Real,2}; #
  minE::Float64 = 1.5e3,
  maxE::Float64 = 0.95 * spec[:BeamEnergy],
  brem::Type{<:NeXLBremsstrahlung} = Castellano2004a,
  mc::Type{<:MatrixCorrection} = Riveros1993,
)

Fit the continuum from ROIs determined from the data within the spectrum (:Composition, :BeamEnergy & :TakeOffAngle). The ROIs are computed using continuumrois(...) and each roi is fit seperately.

source
NeXLSpectrum.detectorresponseFunction
detectorresponse(det::EDSDetector, eff::DetectorEfficiency, incidence::Float64=π/2)::AbstractMatrix

Build a matrix which models the detector response including aspects like the detector efficiency, the resolution, the escape peaks. All the warts that can be modeled within a linear model but not things like coincidence peaks that are non-linear. This function can (!should!) be specialized for more sophisticated detector models that include more warts.

It also dicretizes the input energies on the same scale as the EDSDetector (thus it is square.) This is reasonable when the detector channel width is much less than the resolution.

Example:

genint = computegeneratedintensity(....) # Either characteristic or Bremsstrahlung...
det = simpleEDS(4096, 5.0, 0.0, 132.0, 10)
eff = SDDEfficiency(AP33Tabulation(); thickness=0.0370, deadlayer=30.0e-7, entrance=Film(pure(n"Al"), 10.0e-7))
resp = detectorresponse(det, eff)
# finally compute the measured signal
measured = resp*genint
source
NeXLCore.weightFunction
weight(esc::EscapeArtifact, factor=0.01)

The weight of an EscapeArtifact which is factor * weight(esc.xray).

source
NeXLSpectrum.extentsFunction
extents(cxrs::AbstractVector{<:SpectrumFeature},det::Detector,ampl::Float64)::Vector{UnitRange{Int}}

Determine the contiguous ranges of channels over which the specified collection of X-rays will be measured on the specified detector. The ampl determines the extent of each peak.

source
NeXLSpectrum.profileFunction
profile(energy::Float64, xrayE::Float64, res::Resolution)

Calculates a Gaussian profile for an X-ray of xrayE (eV) for a detector with the specified resolution. Maintains normalization to a sum of unity.

source
profile(ch::Int, xrayE::Float64, det::EDSDetector)

Calculates the profile for an X-ray of xrayE (eV) for a detector with the specified resolution. Performs a crude integration to account for the channel width.

source
NeXLSpectrum.subtractcontinuumFunction
subtractcontinuum(
  spec::Spectrum,
  det::EDSDetector,
  resp::AbstractArray{<:Real,2}; #
  minE::Float64 = 1.5e3,
  maxE::Float64 = 0.95 * spec[:BeamEnergy],
  brem::Type{<:NeXLBremsstrahlung} = Castellano2004a,
  mc::Type{<:MatrixCorrection} = Riveros1993,
)::Spectrum

Computes the characteristic-only spectrum by subtracting the continuum.

source
NeXLSpectrum.heterogeneityFunction
heterogeneity(ffrs::Vector{FilterFitResult}, lbl::ReferenceLabel)

Computes the ratio of the standard deviation of the measured values over the mean calculated uncertainty from the fit. A value near 1 means the sample appears homogeneous and a value greater than 1 means the sample appears heterogeneous.

source

Utility

NeXLSpectrum.drawlineFunction
drawline(func::Function, ci1::CartesianIndex{2}, ci2::CartesianIndex{2}, eachstep::Bool=false)

At each step along the line from ci1 to ci2 call func with a single argument, the row, column coordinates of the pixel as a Tuple{Int, Int}.

If eachstep=true then func is called once when the row changes and once when the column changes. When eachstep=false then func is called less frequently and both row and col may change between calls.

source
NeXLSpectrum.shannon_entropyMethod
shannon_entropy(img::AbstractArray{Gray{N0f8}})

Computes the log2-entropy of the data in the image. The entropy(...) in Images.jl 24.1 is buggy and is removed in 25.0

source
NeXLSpectrum.readSEManticsImageFunction
readSEManticsImage(fn::AbstractString)

Read a SEMantics PNG image and the create an ImageAxes (AxisArray) which associates :x, :y coordinates with the image axes according to the X, Y and FOV in the images.txt file in the same directory.

source
NeXLCore.requiredbutmissingFunction
NeXLCore.requiredbutmissing(ty::Type, spec::Spectrum)

Which properties are missing from spec but are required for the algorithm ty?

source
NeXLCore.hasminrequiredFunction
NeXLCore.hasminrequired(ty::Type, spec::Spectrum)

Does spec have the necessary property items for the algorithm ty?

source
NeXLSpectrum.annotateFunction
annotate(
    img::AxisArray; 
    scale::Bool = true,
    coords::Union{Nothing, AbstractArray{<:AbstractDict{Symbol, <:AbstractFloat}}} = nothing,
    spectra::Union{Nothing, <:AbstractArray{<:Spectrum}} = nothing,
    thumbnails::Union{Nothing, AbstractArray{<:AxisArray}}=nothing,
    labelcoords::Bool = true,
    labelthumbnails::Bool = true,
    marker::Colorant = HSVA(60, 1, 1, 0.3) 
)

Add annotations like scale-bars, acquisition points, image areas to an image and display the result.

  • coords and spectra display as circles (with/without numeric labels).
  • thumbnails display as rectangles overlaying the image
  • scale displays as a scale-bar in the upper-right corner
  • marker is the color of the annotations which is by-default slightly transparent
source
NeXLSpectrum.simulateFunction

" simulate(comp::Material, dose::Float64, e0::Float64, θtoa::Float64, Ω::Float64, det::Detector, resp::Matrix{Float64}; noise=false, vargs...)

Compute a simulated X-ray spectrum for the specified composition material.

Arguments:

  • comp: The Material
  • dose: The electron dose in nA⋅s
  • e0: The beam energy in eV
  • θtoa: The take-off angle in radians
  • Ω: The detector solid angle in steradians
  • det: The detector model
  • resp: The detector response

Returns a Spectrum struct.

source
simulate(
    comp::Material, 
    dose::Float64, 
    e0::Float64, 
    θtoa::Float64, 
    Ω::Float64, 
    det::Detector, 
    resp::Matrix{Float64}; 
    noise=false, 
    vargs...
)

Compute a simulated X-ray spectrum for the specified composition material.

Arguments:

  • comp: The Material
  • dose: The electron dose in nA⋅s
  • e0: The beam energy in eV
  • θtoa: The take-off angle in radians
  • Ω: The detector solid angle in steradians
  • det: The detector model
  • resp: The detector response

Returns a Spectrum struct.

source