609 lines
27 KiB
Matlab
609 lines
27 KiB
Matlab
function DataTable = nappDetermineShiftXps(DataTable, varargin)
|
|
%NAPPDETERMINESHIFTXPS Determine energy shift calibration from XPS images or spectra.
|
|
% This function calculates the energy shift needed to calibrate XPS
|
|
% spectra to a binding energy reference or align all sequence iterations
|
|
% of the same file to a common kinetic energy. For spectra data, it
|
|
% computes absolute binding energy shifts relative to a reference binding
|
|
% energy value. For images, it computes relative kinetic energy shifts
|
|
% between sequence iterations. Shifts are determined by fitting Gaussian
|
|
% peaks or finding maxima in the spectra (simpler option using Matlab's
|
|
% findpeaks function).
|
|
%
|
|
% Needs X-ray energy in data table. Run addPhotonEnergy.m first.
|
|
%
|
|
%
|
|
% SYNTAX
|
|
% ------
|
|
% DataTable = NAPPDETERMINESHIFTXPS(DataTable)
|
|
%
|
|
% DataTable = NAPPDETERMINESHIFTXPS(DataTable, eKinStart, eKinEnd, 'Name', Value)
|
|
%
|
|
% EXAMPLES
|
|
% --------
|
|
% % Determine binding energy shift for C1s spectra using peak finding
|
|
% DataTable = nappDetermineShiftXps(DataTable, ...
|
|
% 'tableVariable', 'spectrum_countsPerSecond', ...
|
|
% 'binEnReference', 285.0, ...
|
|
% 'mode', 'findPeak');
|
|
%
|
|
% % Determine kinetic energy shift for imaging data using the first
|
|
% iteration as reference
|
|
% DataTable = nappDetermineShiftXps(DataTable, 100, 200, ...
|
|
% 'tableVariable', 'image', ...
|
|
% 'kinEnReference', 1, ...
|
|
% 'showFigures', 'on');
|
|
%
|
|
% % Determine binding energy shift using Gaussian fitting with custom
|
|
% parameters specified in tmpFitParameters
|
|
% FitPar.fitMode = 'position';
|
|
% FitPar.refPeak = 1;
|
|
% DataTable = nappDetermineShiftXps(DataTable, 280, 290, ...
|
|
% 'tableVariable', 'spectrum_countsPerSecond', ...
|
|
% 'binEnReference', 284.8, ...
|
|
% 'mode', 'gaussFit', ...
|
|
% 'FitParameter', tmpFitParameters);
|
|
%
|
|
% tmpFitParameters.peakArea = [1e5 0.5 0.1 0.1 0.1] ;
|
|
% tmpFitParameters.peakWidth = 1.2;
|
|
% tmpFitParameters.peakPos = [198 1.6 10 10 10];
|
|
%
|
|
% tmpFitParameters.peakAreaUpperBond = [Inf 0.5 0.1 0.1 0.1];
|
|
% tmpFitParameters.peakWidthUpperBond = [3 NaN NaN NaN NaN];
|
|
% tmpFitParameters.peakPosUpperBond = [200 1.6 10 10 10];
|
|
%
|
|
% tmpFitParameters.peakAreaLowerBond = [100 10 0.5 0.1 0.1 0.1];
|
|
% tmpFitParameters.peakWidthLowerBond = [1.2 NaN NaN NaN NaN];
|
|
% tmpFitParameters.peakPosLowerBxxond = [195 1.3 0.1 0.1 0.1];
|
|
%
|
|
% tmpFitParameters.backgroundRange = [6 5];
|
|
% tmpFitParameters.weight = [198 200 0.1];
|
|
% tmpFitParameters.fitMode = 'widthPosition';
|
|
%
|
|
% PROCESSING OVERVIEW
|
|
% -------------------
|
|
% The function operates differently depending on data type:
|
|
%
|
|
% FOR IMAGES (tableVariable = 'image' or 'imageRoh'):
|
|
% 1. For each image in tableSelection:
|
|
% - Extract kinetic energy axis (scientaEkin_eV or scientaEkinRoh_eV)
|
|
% - Fit 2-peak Gaussian model to each sequence iteration (WIP: change
|
|
% this to use nappFitGaussXPS.m)
|
|
% - Extract peak position from fitting results
|
|
% - Normalize shifts relative to reference iteration (kinEnReference)
|
|
% 2. Result: kineticEnergyShift column contains shifts relative to reference iteration
|
|
% (reference iteration has shift = 0)
|
|
%
|
|
% FOR SPECTRA (tableVariable = 'spectrum_counts', 'spectrum_countsPerSecond', or 'spectrum_countsNorm'):
|
|
% 1. For each spectrum in tableSelection:
|
|
% - Optional: restrict energy range using eKinStart/eKinEnd
|
|
% - Find peak position using either:
|
|
% * findPeak method: standard matlab peak detection
|
|
% * gaussFit method: Gaussian fitting with global optimization
|
|
% - Compute binding energy shift as: xRayPhotonEnergy - peakPosition - binEnReference
|
|
% 2. Result: bindingEnergyShift column contains absolute shifts relative to reference BE
|
|
% (positive shift means peak is at lower kinetic energy than expected)
|
|
%
|
|
% PEAK DETECTION METHODS
|
|
% ----------------------
|
|
% Two methods are available for finding peak positions in spectral
|
|
%
|
|
% 1. 'findPeak' (default): Fast automated detection using MATLAB's findpeaks
|
|
% function with moving mean smoothing. Uses minPeakWidth, minPeakHeight,
|
|
% and movMeanWindow parameters to control sensitivity. Best for well-defined
|
|
% peaks with good signal-to-noise ratio. Accuracy: ±0.1-0.2 eV typical.
|
|
%
|
|
% 2. 'gaussFit': Slower but more accurate Gaussian fitting with global
|
|
% optimization using GlobalSearch. Fits background, Gaussian peaks,
|
|
% and extracts precise peak positions. Requires FitParameter structure.
|
|
% Best for noisy data or when high precision is required. Accuracy: ±0.05 eV typical.
|
|
%
|
|
% IMPORTANT LIMITATION: This function requires distinct photoelectron peaks
|
|
% to determine energy shifts. It will NOT work for background-only spectra
|
|
% or featureless data. The fitting algorithms assume at least one measurable
|
|
% peak is present in the specified energy range.
|
|
%
|
|
% INPUT ARGUMENTS
|
|
% ---------------
|
|
% DataTable table
|
|
% Data table containing XPS image or spectral data. Must contain:
|
|
% - scientaEkin_eV or scientaEkinRoh_eV: kinetic energy axis
|
|
% - tableVariable data: images or spectra
|
|
% - xRayEkin_eV: X-ray photon energy (eV) for each row
|
|
%
|
|
% eKinStart numeric scalar (optional)
|
|
% Lower kinetic energy bound for peak fitting (eV). If omitted, uses
|
|
% first energy point. Restricts the fitting region to improve accuracy.
|
|
%
|
|
% eKinEnd numeric scalar (optional)
|
|
% Upper kinetic energy bound for peak fitting (eV). If omitted, uses
|
|
% last energy point. Restricts the fitting region to improve accuracy.
|
|
%
|
|
% NAME-VALUE PAIR ARGUMENTS (optional)
|
|
% ------------------------------------
|
|
% 'tableSelection' logical array (default: all rows)
|
|
% Logical index to filter which rows in DataTable to process.
|
|
% Example: contains(DataTable.name, 'C1s') to process only C1s spectra
|
|
%
|
|
% 'tableVariable' char vector (default: 'spectrum_countsPerSecond')
|
|
% Data type to process. Must be one of:
|
|
%
|
|
% Spectral data (compute binding energy shifts):
|
|
% - 'spectrum_counts'
|
|
% - 'spectrum_countsPerSecond'
|
|
% - 'spectrum_countsNorm'
|
|
%
|
|
% Imaging data (compute kinetic energy shifts):
|
|
% - 'image'
|
|
% - 'imageRoh'
|
|
%
|
|
% 'mode' char vector (default: 'findPeak')
|
|
% Peak detection method (spectral data only):
|
|
% - 'findPeak': automated peak detection using moving mean smoothing
|
|
% - 'gaussFit': 2-peak Gaussian fitting with global optimization
|
|
%
|
|
% 'binEnReference' numeric scalar (required for spectra)
|
|
% Reference binding energy value (eV) for spectral calibration.
|
|
% Typically the known BE of the primary peak (e.g., 285.0 eV for C1s C-H).
|
|
% Example: 'binEnReference', 285.0
|
|
%
|
|
% 'kinEnReference' numeric scalar (default: 1)
|
|
% Reference iteration index for imaging data. All kinetic energy shifts
|
|
% are computed relative to this iteration. Reference iteration gets shift = 0.
|
|
% Range: 1 to number of image columns.
|
|
%
|
|
% 'mode' char vector (default: 'findPeak')
|
|
% Peak finding method for spectra:
|
|
% - 'findPeak': uses findpeaks with smoothing and height/width constraints
|
|
% - 'gaussFit': uses nappFitGaussXps with global optimization (slower, more accurate)
|
|
%
|
|
% 'minPeakWidth' numeric scalar (default: 0.1)
|
|
% Minimum peak width (eV) for findPeak method. Peaks narrower than this
|
|
% are rejected. Typical range: 0.05-0.5 eV depending on spectrometer resolution.
|
|
%
|
|
% 'minPeakHeight' numeric scalar (default: 0.1)
|
|
% Minimum peak height (intensity) for findPeak method. Peaks below this
|
|
% are rejected. Units match tableVariable data.
|
|
%
|
|
% 'movMeanWindow' numeric scalar (default: 1)
|
|
% Moving mean window width (number of points) for smoothing before peak
|
|
% detection in findPeak method. Reduces noise; typical range: 1-5 points.
|
|
%
|
|
% 'FitParameter' struct
|
|
% Structure containing Gaussian fitting parameters for gaussFit mode.
|
|
% Required fields:
|
|
% - fitMode: 'position', 'width', 'unconstrained', etc.
|
|
% - refPeak: peak number (1 or 2) whose position is used for shift calculation
|
|
% Example: FitPar.fitMode = 'position'; FitPar.refPeak = 1;
|
|
%
|
|
% 'refPeak' numeric scalar (default: 1)
|
|
% Peak number (1 or 2) whose position is used for binding energy shift
|
|
% calculation in gaussFit mode. Peak 1 is typically the main photoelectron peak.
|
|
%
|
|
% 'showFigures' char vector (default: 'off')
|
|
% Display fitting results:
|
|
% - 'on': show Gaussian fit plot, peak table, and peak finding results
|
|
% - 'off': suppress all figure output
|
|
%
|
|
% OUTPUT ARGUMENTS
|
|
% ----------------
|
|
% DataTable table
|
|
% Updated data table with new columns:
|
|
%
|
|
% For imaging data (tableVariable = 'image' or 'imageRoh'):
|
|
% - kineticEnergyShift: column vector of relative kinetic energy shifts (eV)
|
|
% relative to kinEnReference iteration. Reference = 0.
|
|
%
|
|
% For spectral data (tableVariable = 'spectrum_*'):
|
|
% - bindingEnergyShift: column vector of absolute binding energy shifts (eV)
|
|
% relative to binEnReference. Typically small values (±0.5 eV range).
|
|
%
|
|
% For gaussFit mode on spectra, additional columns are added:
|
|
% - peakPosition: [peak1, peak2] positions (eV) for each spectrum
|
|
% - peakArea: [peak1, peak2] areas (eV·counts) for each spectrum
|
|
% - peakFWHM: [peak1, peak2] full-width half-maximum (eV) for each spectrum
|
|
% - peakBackground: linear background parameters [slope, intercept]
|
|
%
|
|
% NOTES
|
|
% -----
|
|
% - PREREQUISITE: X-ray photon energy (xRayEkin_eV) must be present in DataTable.
|
|
% If empty, function returns with warning.
|
|
% - For imaging: Gaussian fitting assumes 2 peaks in each iteration. Adjust
|
|
% eKinStart/eKinEnd to isolate peak region if data is complex.
|
|
% - For spectra: findPeak method is fast but less accurate (±0.1-0.2 eV typical).
|
|
% gaussFit method is slower but more accurate (±0.05 eV typical).
|
|
% - Peak detection may fail if signal-to-noise is poor. Check showFigures output
|
|
% to diagnose fitting quality.
|
|
% - For imageRoh variant: uses scientaEkinRoh_eV if available; falls back to
|
|
% scientaEkin_eV if not populated.
|
|
% - Binding energy shifts can be negative (peak at higher BE than reference) or
|
|
% positive (peak at lower BE than reference).
|
|
% - Figure windows are automatically closed after processing to prevent clutter.
|
|
%
|
|
% SEE ALSO
|
|
% --------
|
|
% nappCopyShift, nappBindingEnergyXps, nappFitGaussXps, nappAlignImages,
|
|
% findpeaks, GaussFitOnePeak
|
|
defaultKinEnReference = 1;
|
|
defaultTableVariable = 'spectrum_countsPerSecond';
|
|
expectedTableVariable = {'spectrum_counts','spectrum_countsPerSecond','spectrum_countsNorm','image','imageRoh'};
|
|
defaultTableSelection = true(size(DataTable,1));
|
|
defaultMinPeakWidth = 0.1;
|
|
defaultMinPeakHeight = 0.1;
|
|
defaultMovMeanWindow = 1;
|
|
defaultRefPeak = 1;
|
|
defaultShowFigures = 'off';
|
|
expectedShowFigures = {'on','off'};
|
|
defaultMode = 'findPeak';
|
|
expectedMode = {'gaussFit','findPeak'};
|
|
% defaultFitMode = 'widthPosition';
|
|
% expectedFitMode = {'unconstrained','widthPosition','position','width','unconstrainedAndWidthPosition'};
|
|
|
|
|
|
%%
|
|
|
|
% Parse input
|
|
p = inputParser;
|
|
addRequired(p,'DataTable',@istable);
|
|
addOptional(p,'eKinStart',@isnumeric)
|
|
addOptional(p,'eKinEnd',@isnumeric)
|
|
addParameter(p,'FitParameter',@isstruct);
|
|
addParameter(p,'binEnReference',@isnumeric)
|
|
addParameter(p,'kinEnReference',defaultKinEnReference,@isnumeric)
|
|
addParameter(p,'tableSelection',defaultTableSelection,@islogical);
|
|
addParameter(p,'minPeakWidth',defaultMinPeakWidth,@isnumeric);
|
|
addParameter(p,'minPeakHeight',defaultMinPeakHeight,@isnumeric);
|
|
addParameter(p,'movMeanWindow',defaultMovMeanWindow ,@isnumeric);
|
|
addParameter(p,'tableVariable',defaultTableVariable,...
|
|
@(x) any(validatestring(x,expectedTableVariable)));
|
|
addParameter(p,'showFigures',defaultShowFigures,...
|
|
@(x) any(validatestring(x,expectedShowFigures)));
|
|
addParameter(p,'mode',defaultMode, @(x) any(validatestring(x,expectedMode)));
|
|
% addParameter(p,'fitMode',defaultFitMode, @(x) any(validatestring(x,expectedFitMode)));
|
|
addParameter(p,'refPeak',defaultRefPeak,@isnumeric);
|
|
parse(p,DataTable,varargin{:});
|
|
|
|
|
|
% check input (all elements of input variable need to have non-zero size)
|
|
if ~all(cell2mat(cellfun(@size,p.Results.DataTable(p.Results.tableSelection,:).(p.Results.tableVariable),'UniformOutput',false)),'all')
|
|
disp('The table variable appears to be empty. I quit now')
|
|
error('The table variable appears to be empty.')
|
|
end
|
|
|
|
% add empty variable to table
|
|
if contains(p.Results.tableVariable,'image')
|
|
existColumn = strcmp('kineticEnergyShift',DataTable.Properties.VariableNames);
|
|
if ~any(existColumn)
|
|
DataTable.kineticEnergyShift = cell(size(DataTable,1),1);
|
|
end
|
|
else
|
|
existColumn = strcmp('bindingEnergyShift',DataTable.Properties.VariableNames);
|
|
if ~any(existColumn)
|
|
DataTable.bindingEnergyShift = NaN(size(DataTable,1),1);
|
|
end
|
|
end
|
|
%Colors; define PSI cooperate identity colours
|
|
PsiFarbe.hellgrau = [226 226 226]/255;
|
|
PsiFarbe.mittelgrau = [184 184 184]/255;
|
|
PsiFarbe.dunkelgrau = [78 78 77]/255;
|
|
PsiFarbe.blaugrau = [236 236 243]/255; % This is actually not a PSI colour
|
|
PsiFarbe.gelb = [239 200 37]/255;
|
|
PsiFarbe.hellorange = [227 164 31]/255;
|
|
PsiFarbe.dunkelorange = [207 92 22]/255;
|
|
PsiFarbe.hellrot = [164 36 23]/255;
|
|
PsiFarbe.dunkelbraun = [114 85 65]/255;
|
|
PsiFarbe.hellbraun = [127 119 47]/255;
|
|
PsiFarbe.mittelgruen = [135 150 59]/255;
|
|
PsiFarbe.hellgruen = [195 214 47]/255;
|
|
PsiFarbe.dunkelgruen = [72 130 33]/255;
|
|
PsiFarbe.dunkelrot = [104 47 86]/255;
|
|
PsiFarbe.dunkelblau = [37 65 108]/255;
|
|
PsiFarbe.mittelgruenblau = [61 164 168]/255;
|
|
PsiFarbe.hellblau = [58 162 190]/255;
|
|
PsiFarbe.mittelblau = [50 118 158]/255;
|
|
|
|
% confirm that X-ray energy is set
|
|
if any(isnan(cell2mat(p.Results.DataTable(p.Results.tableSelection,:).xRayEkin_eV)))
|
|
disp('Opps, X-ray energy is empty. Please add manually. Will stop now');
|
|
return
|
|
end
|
|
if matches(p.Results.tableVariable,'image')
|
|
% repeat this for all columns in input.
|
|
for ii = 1:size(p.Results.DataTable(p.Results.tableSelection,:),1)
|
|
%%
|
|
try
|
|
if isnumeric(p.Results.eKinStart)
|
|
start = find(p.Results.DataTable(p.Results.tableSelection,:).scientaEkin_eV{ii}...
|
|
>= p.Results.eKinStart,1);
|
|
stop = find(p.Results.DataTable(p.Results.tableSelection,:).scientaEkin_eV{ii}...
|
|
>= p.Results.eKinEnd,1);
|
|
else
|
|
start = 1;
|
|
stop = length(p.Results.DataTable(p.Results.tableSelection,:).scientaEkin_eV{ii});
|
|
end
|
|
catch
|
|
warning('There is something wrong with the eKinStart and/or eKinEnd parameters');
|
|
disp('There is something wrong with the eKinStart and/or eKinEnd parameters');
|
|
end
|
|
|
|
for ij = 1:size(p.Results.DataTable(p.Results.tableSelection,:).(p.Results.tableVariable){ii},2)
|
|
[B,x1,y1] = GaussFitOnePeak([p.Results.DataTable(p.Results.tableSelection,:).scientaEkin_eV{ii}(start:stop)',...
|
|
p.Results.DataTable(p.Results.tableSelection,:).(p.Results.tableVariable){ii}(start:stop,ij)]);
|
|
tmp_1(ij) = p.Results.DataTable(p.Results.tableSelection,:).xRayEkin_eV{ii}-B(3);
|
|
% plot result
|
|
if matches(p.Results.showFigures,'on')
|
|
% define equations to be plotted
|
|
Gauss = @(x,Hgt, Pos, Wdt) Hgt./(Wdt.*sqrt(pi/2))*exp(-2.*((x - Pos).^2/Wdt.^2));
|
|
linear = @(x,slp,intrcpt) slp*x+intrcpt;
|
|
figure( 'Name', 'Gauss Multi Fit ')
|
|
tiledlayout(1,1)
|
|
nexttile;
|
|
%Define gauss peaks with hight, position, width from fitting results
|
|
yData_pk1= Gauss(x1 ,B(1), B(3), B(2));
|
|
yData_pk2= Gauss(x1 ,B(4), B(6), B(5));
|
|
bkgrd = linear(x1,B(7),B(8));
|
|
% plot data
|
|
plot(x1,y1,'.')
|
|
hold on
|
|
% plot fit
|
|
h = plot(x1,yData_pk1+yData_pk2+bkgrd);
|
|
h.LineWidth=2;
|
|
h.Color = PsiFarbe.dunkelblau;
|
|
% plot each gauss peak
|
|
h = plot(x1,yData_pk1+bkgrd);
|
|
h.LineWidth=1.5;
|
|
h.Color = PsiFarbe.gelb;
|
|
h = plot(x1,yData_pk2+bkgrd);
|
|
h.LineWidth=1.5;
|
|
h.Color = PsiFarbe.gelb;
|
|
|
|
legend OFF
|
|
xlabel('Kinetic Energy [eV]');
|
|
ylabel('PE intensity [-]');
|
|
grid on
|
|
|
|
% List results in table
|
|
|
|
disp('Results');
|
|
disp(' ' );
|
|
disp(cell2table({sprintf('%.1f',B(3)),...
|
|
sprintf('%.1f',B(2)*sqrt(log(4))),...
|
|
sprintf('%.1f',B(2)),...
|
|
sprintf('%7.0f',B(1))},...
|
|
'VariableNames',{'position';'FWHW';'width';'area'},...
|
|
'RowNames',{'pk 1'}));
|
|
end
|
|
end
|
|
%%
|
|
|
|
|
|
DataTable(p.Results.tableSelection,:).kineticEnergyShift{ii} = tmp_1-tmp_1(p.Results.kinEnReference);
|
|
clearvars tmp_1
|
|
fh = findobj( 'Type', 'Figure', 'Name', 'Gauss Multi Fit ');
|
|
close(fh)
|
|
fh = findobj( 'Type', 'Figure', 'Name', 'GlobalSearch');
|
|
close(fh)
|
|
end
|
|
elseif matches(p.Results.tableVariable,'imageRoh')
|
|
% repeat this for all columns in input.
|
|
for ii = 1:size(p.Results.DataTable(p.Results.tableSelection,:),1)
|
|
|
|
%%
|
|
try
|
|
if isnumeric(p.Results.eKinStart) && isempty(DataTable(p.Results.tableSelection,:).scientaEkinRoh_eV{ii})
|
|
start = find(p.Results.DataTable(p.Results.tableSelection,:).scientaEkin_eV{ii}...
|
|
>= p.Results.eKinStart,1);
|
|
stop = find(p.Results.DataTable(p.Results.tableSelection,:).scientaEkin_eV{ii}...
|
|
>= p.Results.eKinEnd,1);
|
|
elseif isnumeric(p.Results.eKinStart) && ~isempty(DataTable(p.Results.tableSelection,:).scientaEkinRoh_eV{ii})
|
|
start = find(p.Results.DataTable(p.Results.tableSelection,:).scientaEkinRoh_eV{ii}...
|
|
>= p.Results.eKinStart,1);
|
|
stop = find(p.Results.DataTable(p.Results.tableSelection,:).scientaEkinRoh_eV{ii}...
|
|
>= p.Results.eKinEnd,1);
|
|
elseif isempty(DataTable(p.Results.tableSelection,:).scientaEkinRoh_eV{ii})
|
|
start = 1;
|
|
stop = length(p.Results.DataTable(p.Results.tableSelection,:).scientaEkin_eV{ii});
|
|
else
|
|
start = 1;
|
|
stop = length(DataTable(p.Results.tableSelection,:).scientaEkinRoh_eV{ii});
|
|
end
|
|
catch
|
|
warning('There is something wrong with the eKinStart and/or eKinEnd parameters');
|
|
disp('There is something wrong with the eKinStart and/or eKinEnd parameters');
|
|
end
|
|
for ij = 1:size(p.Results.DataTable(p.Results.tableSelection,:).(p.Results.tableVariable){ii},2)
|
|
if isempty(DataTable(p.Results.tableSelection,:).scientaEkinRoh_eV{ii})
|
|
[B,x1,y1] = GaussFitOnePeak([p.Results.DataTable(p.Results.tableSelection,:).scientaEkin_eV{ii}(start:stop)',...
|
|
p.Results.DataTable(p.Results.tableSelection,:).(p.Results.tableVariable){ii}(start:stop,ij)]);
|
|
else
|
|
[B,x1,y1] = GaussFitOnePeak([DataTable(p.Results.tableSelection,:).scientaEkinRoh_eV{ii}(start:stop)',...
|
|
p.Results.DataTable(p.Results.tableSelection,:).(p.Results.tableVariable){ii}(start:stop,ij)]);
|
|
end
|
|
|
|
tmp_1(ij) = p.Results.DataTable(p.Results.tableSelection,:).xRayEkin_eV{ii}-B(3);
|
|
% plot result
|
|
if matches(p.Results.showFigures,'on')
|
|
% define equations to be plotted
|
|
Gauss = @(x,Hgt, Pos, Wdt) Hgt./(Wdt.*sqrt(pi/2))*exp(-2.*((x - Pos).^2/Wdt.^2));
|
|
linear = @(x,slp,intrcpt) slp*x+intrcpt;
|
|
figure( 'Name', 'Gauss Multi Fit ')
|
|
tiledlayout(1,1)
|
|
nexttile;
|
|
%Define gauss peaks with hight, position, width from fitting results
|
|
yData_pk1= Gauss(x1 ,B(1), B(3), B(2));
|
|
yData_pk2= Gauss(x1 ,B(4), B(6), B(5));
|
|
bkgrd = linear(x1,B(7),B(8));
|
|
% plot data
|
|
plot(x1,y1,'.')
|
|
hold on
|
|
% plot fit
|
|
h = plot(x1,yData_pk1+yData_pk2+bkgrd);
|
|
h.LineWidth=2;
|
|
h.Color = PsiFarbe.dunkelblau;
|
|
% plot each gauss peak
|
|
h = plot(x1,yData_pk1+bkgrd);
|
|
h.LineWidth=1.5;
|
|
h.Color = PsiFarbe.gelb;
|
|
h = plot(x1,yData_pk2+bkgrd);
|
|
h.LineWidth=1.5;
|
|
h.Color = PsiFarbe.gelb;
|
|
|
|
legend OFF
|
|
xlabel('Kinetic Energy [eV]');
|
|
ylabel('PE intensity [-]');
|
|
grid on
|
|
|
|
% List results in table
|
|
|
|
disp('Results');
|
|
disp(' ' );
|
|
disp(cell2table({sprintf('%.1f',B(3)),...
|
|
sprintf('%.1f',B(2)*sqrt(log(4))),...
|
|
sprintf('%.1f',B(2)),...
|
|
sprintf('%7.0f',B(1))},...
|
|
'VariableNames',{'position';'FWHW';'width';'area'},...
|
|
'RowNames',{'pk 1'}));
|
|
end
|
|
end
|
|
%%
|
|
|
|
DataTable(p.Results.tableSelection,:).kineticEnergyShift{ii} = tmp_1-tmp_1(p.Results.kinEnReference);
|
|
clearvars tmp_1
|
|
fh = findobj( 'Type', 'Figure', 'Name', 'Gauss Multi Fit ');
|
|
close(fh)
|
|
fh = findobj( 'Type', 'Figure', 'Name', 'GlobalSearch');
|
|
close(fh)
|
|
end
|
|
else
|
|
% repeat this for all columns in input.
|
|
if matches(p.Results.mode,'findPeak')
|
|
for ii = 1:size(p.Results.DataTable(p.Results.tableSelection,:),1)
|
|
|
|
%%
|
|
try
|
|
if isnumeric(p.Results.eKinStart)
|
|
start = find(p.Results.DataTable(p.Results.tableSelection,:).scientaEkin_eV{ii}...
|
|
>= p.Results.eKinStart,1);
|
|
stop = find(p.Results.DataTable(p.Results.tableSelection,:).scientaEkin_eV{ii}...
|
|
>= p.Results.eKinEnd,1);
|
|
else
|
|
start = 1;
|
|
stop = length(p.Results.DataTable(p.Results.tableSelection,:).scientaEkin_eV{ii});
|
|
end
|
|
catch
|
|
warning('There is something wrong with the eKinStart and/or eKinEnd parameters');
|
|
disp('There is something wrong with the eKinStart and/or eKinEnd parameters');
|
|
end
|
|
|
|
if contains(p.Results.showFigures,'on')
|
|
figure
|
|
|
|
findpeaks(movmean(p.Results.DataTable(p.Results.tableSelection,:).(p.Results.tableVariable){ii}(start:stop),...
|
|
p.Results.movMeanWindow),...
|
|
p.Results.DataTable(p.Results.tableSelection,:).scientaEkin_eV{ii}(start:stop),...
|
|
'MinPeakWidth',p.Results.minPeakWidth,'MinPeakHeight',p.Results.minPeakHeight)
|
|
end
|
|
[pks,locs]=findpeaks(movmean(p.Results.DataTable(p.Results.tableSelection,:).(p.Results.tableVariable){ii}(start:stop),...
|
|
p.Results.movMeanWindow),...
|
|
p.Results.DataTable(p.Results.tableSelection,:).scientaEkin_eV{ii}(start:stop),...
|
|
'MinPeakWidth',p.Results.minPeakWidth,'MinPeakHeight',p.Results.minPeakHeight);
|
|
|
|
%%
|
|
|
|
[~,I] = max(pks);
|
|
tmp_1(ii,1) = p.Results.DataTable(p.Results.tableSelection,:).xRayEkin_eV{ii}-locs(I)-p.Results.binEnReference;
|
|
clearvars tmp_2
|
|
end
|
|
DataTable.bindingEnergyShift(p.Results.tableSelection) = tmp_1;
|
|
|
|
elseif matches(p.Results.mode,'gaussFit')
|
|
tmpFitResults = nappFitGaussXps(DataTable,p.Results.FitParameter,p.Results.tableSelection,...
|
|
'energyScale','kinetic','mode',p.Results.FitParameter.fitMode,'showFigures','on','tableVariable',p.Results.tableVariable);
|
|
|
|
DataTable.bindingEnergyShift(p.Results.tableSelection) = tmpFitResults(p.Results.tableSelection,:).xRayEkinRange_eV(:,1) - tmpFitResults(p.Results.tableSelection,:).peakPosition(:,p.Results.FitParameter.refPeak) - repmat(p.Results.binEnReference,size(DataTable.bindingEnergyShift(p.Results.tableSelection)));
|
|
DataTable.peakPosition(p.Results.tableSelection,:) = tmpFitResults(p.Results.tableSelection,:).peakPosition;
|
|
DataTable.peakArea(p.Results.tableSelection,:) = tmpFitResults(p.Results.tableSelection,:).peakArea;
|
|
DataTable.peakFWHM(p.Results.tableSelection,:) = tmpFitResults(p.Results.tableSelection,:).peakFWHM;
|
|
% added 6.2025 tbr
|
|
DataTable.peakBackground(p.Results.tableSelection,:) = tmpFitResults(p.Results.tableSelection,:).peakBackground;
|
|
%%
|
|
|
|
|
|
|
|
end
|
|
|
|
|
|
end
|
|
try
|
|
fh = findobj( 'Type', 'Figure', 'Name', 'Gauss Multi Fit ');
|
|
close(fh)
|
|
catch
|
|
end
|
|
try
|
|
fh = findobj( 'Type', 'Figure', 'Name', 'GlobalSearch');
|
|
close(fh)
|
|
catch
|
|
end
|
|
|
|
end
|
|
function [B,x1,y1] = GaussFitOnePeak(data1)
|
|
|
|
|
|
[x1, y1] = prepareCurveData(data1(:,1),data1(:,2));
|
|
|
|
% make data slimmer so that background makes sense
|
|
[~,I] = max(y1);
|
|
tmp1 = find(x1>=x1(I)-8,1);
|
|
tmp2 = find(x1>=x1(I)+5,1);
|
|
if isempty(tmp2)
|
|
y1 = y1(tmp1:end);
|
|
x1 = x1(tmp1:end);
|
|
else
|
|
y1 = y1(tmp1:tmp2);
|
|
x1 = x1(tmp1:tmp2);
|
|
end
|
|
[~,I] = max(y1);
|
|
|
|
% Define equation to fit data
|
|
% I use the area Gauss equations here, which I took from the OriginPro
|
|
% documentation. https://www.originlab.com/doc/Origin-Help/Gauss-FitFunc
|
|
%
|
|
% Gauss = area/(width*sqrt(pi/2)) * exp(-2*(x-position)^2/width^2)
|
|
%
|
|
% where FWHM = width * sqrt(ln(4))
|
|
|
|
% Data set 1:
|
|
yfcn1 = @(b,x)...
|
|
b(1)/(b(2)*sqrt(pi/2))*exp(-2.*((x(:,1) - b(3)).^2/b(2).^2)) ... % Gauss peak 1: C-H, Ac-
|
|
+ b(4)/(b(5)*sqrt(pi/2))*exp(-2.*((x(:,1)- b(6)).^2/b(5).^2))... % Gauss peak 2: C-OOH, Ac-
|
|
+b(7).*x(:,1)+b(8); % linear background
|
|
|
|
B0 = [1e5 1.5 x1(I) ...
|
|
1e5 1.0 x1(I)-4 ...
|
|
1 0];
|
|
|
|
|
|
% lower bonds
|
|
lb = [100 0.5 x1(I)-1. ...
|
|
100 1.0 x1(I)-5 ...
|
|
-Inf -Inf];
|
|
|
|
|
|
% upper bond
|
|
ub = [Inf 3 x1(I)+1 ...
|
|
Inf 5 x1(I)-1 ...
|
|
Inf Inf];
|
|
|
|
|
|
% Fit equations: Gauss eq. 1 with data
|
|
yfcn = @(b,x) yfcn1(b,x);
|
|
% We find an optimal fit by minimising the residuals:
|
|
fitfcnw = @(b) sum((y1 - yfcn(b,x1)).^2);
|
|
|
|
problem = createOptimProblem('fmincon', 'x0',B0, 'lb',lb,'ub',ub, 'objective',fitfcnw);
|
|
gs = GlobalSearch('PlotFcns',@gsplotbestf);
|
|
B = run(gs,problem);
|
|
|
|
|
|
|
|
|
|
end |