How to calculate the width where the intensity is higher than 20% of the maximum intensity?

3 ビュー (過去 30 日間)
To do this, I think I have to utilize the findpeaks function. Then interpolate the 20% to the full width at half maximum (FWHM). So we can find its width?
Below is the code for an attempt to get the max intensity and its FWHM:
[pks,widths] = findpeaks(Int,x2,'NPeaks',1,'WidthReference','halfheight');

採用された回答

Mathieu NOE
Mathieu NOE 2022 年 4 月 29 日
hello Angelo
I have no doubt you can use findpeaks
but sometimes I prefer my own way.
see demo below - this is to find W20 and extract the portion of data between the two crossing points
NB that the function crossing_V7 use linear interpolation which is more accurate (on the x axis) then simply finding the first x data element above 20% (this depends of your x spacing)
hope it helps
clc
clearvars
% data
load('position.mat');
x = x2; clear x2
load('Intensity.mat');
y = Int; clear Int
threshold = 0.2*max(y); % 20% of peak amplitude
[t0_pos1,s0_pos1,t0_neg1,s0_neg1]= crossing_V7(y,x,threshold,'linear'); % positive (pos) and negative (neg) slope crossing points
% ind => time index (samples)
% t0 => corresponding time (x) values
% s0 => corresponding function (y) values , obviously they must be equal to "threshold"
% depending of the threshold level we have either one or three peaks
% if we have 3 peaks then pick the second one
if numel(s0_pos1) >1
s0_pos1 = s0_pos1(2);
s0_neg1 = s0_neg1(2);
t0_pos1 = t0_pos1(2);
t0_neg1 = t0_neg1(2);
end
% x difference between t0_neg and t0_pos = W20
W20 = t0_neg1 - t0_pos1
% extract the data and do interpolation
xx = linspace(t0_pos1,t0_neg1,100);
yy = interp1(x,y,xx);
figure(1)
plot(x,y,'b',t0_pos1,s0_pos1,'dr',t0_neg1,s0_neg1,'dg',xx,yy,'linewidth',2,'markersize',12);grid on
legend('signal','signal positive slope crossing points','signal negative slope crossing points','selected data');
xlim([-6e-4 6e-4]);
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function [t0_pos,s0_pos,t0_neg,s0_neg] = crossing_V7(S,t,level,imeth)
% [ind,t0,s0,t0close,s0close] = crossing_V6(S,t,level,imeth,slope_sign) % older format
% CROSSING find the crossings of a given level of a signal
% ind = CROSSING(S) returns an index vector ind, the signal
% S crosses zero at ind or at between ind and ind+1
% [ind,t0] = CROSSING(S,t) additionally returns a time
% vector t0 of the zero crossings of the signal S. The crossing
% times are linearly interpolated between the given times t
% [ind,t0] = CROSSING(S,t,level) returns the crossings of the
% given level instead of the zero crossings
% ind = CROSSING(S,[],level) as above but without time interpolation
% [ind,t0] = CROSSING(S,t,level,par) allows additional parameters
% par = {'none'|'linear'}.
% With interpolation turned off (par = 'none') this function always
% returns the value left of the zero (the data point thats nearest
% to the zero AND smaller than the zero crossing).
%
% check the number of input arguments
error(nargchk(1,4,nargin));
% check the time vector input for consistency
if nargin < 2 | isempty(t)
% if no time vector is given, use the index vector as time
t = 1:length(S);
elseif length(t) ~= length(S)
% if S and t are not of the same length, throw an error
error('t and S must be of identical length!');
end
% check the level input
if nargin < 3
% set standard value 0, if level is not given
level = 0;
end
% check interpolation method input
if nargin < 4
imeth = 'linear';
end
% make row vectors
t = t(:)';
S = S(:)';
% always search for zeros. So if we want the crossing of
% any other threshold value "level", we subtract it from
% the values and search for zeros.
S = S - level;
% first look for exact zeros
ind0 = find( S == 0 );
% then look for zero crossings between data points
S1 = S(1:end-1) .* S(2:end);
ind1 = find( S1 < 0 );
% bring exact zeros and "in-between" zeros together
ind = sort([ind0 ind1]);
% and pick the associated time values
t0 = t(ind);
s0 = S(ind);
if ~isempty(ind)
if strcmp(imeth,'linear')
% linear interpolation of crossing
for ii=1:length(t0)
%if abs(S(ind(ii))) >= eps(S(ind(ii))) % MATLAB V7 et +
if abs(S(ind(ii))) >= eps*abs(S(ind(ii))) % MATLAB V6 et - EPS * ABS(X)
% interpolate only when data point is not already zero
NUM = (t(ind(ii)+1) - t(ind(ii)));
DEN = (S(ind(ii)+1) - S(ind(ii)));
slope = NUM / DEN;
slope_sign(ii) = sign(slope);
t0(ii) = t0(ii) - S(ind(ii)) * slope;
s0(ii) = level;
end
end
end
% extract the positive slope crossing points
ind_pos = find(sign(slope_sign)>0);
t0_pos = t0(ind_pos);
s0_pos = s0(ind_pos);
% extract the negative slope crossing points
ind_neg = find(sign(slope_sign)<0);
t0_neg = t0(ind_neg);
s0_neg = s0(ind_neg);
else
% empty output
ind_pos = [];
t0_pos = [];
s0_pos = [];
% extract the negative slope crossing points
ind_neg = [];
t0_neg = [];
s0_neg = [];
end
end
  2 件のコメント
Lope
Lope 2022 年 5 月 2 日
@Mathieu NOE Thank you so much! This really helps. Thank you also for the extra effort in giving explanations for the entire process.

サインインしてコメントする。

その他の回答 (0 件)

カテゴリ

Help Center および File ExchangeParticle & Nuclear Physics についてさらに検索

タグ

製品


リリース

R2021a

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by