Peak area through Findpeak

Hello Everyone,

I have the following data: as attched m file. I can fine the peak location, width through peakfind command. Now i want to find the area under the peak. would you please guide me how i can find using findpeaks.

[pks20_4,loc20_4,w20_4,p20_4] =findpeaks(Y20_4,X4);

Thank you,

Answer 1

[Star Strider]

It depends on how you define ‘peak’.

LD1 = load('X4.mat');
X4 = LD1.X4;
LD2 = load('Y20_4.mat');
Y20_4 = LD2.Y20_4;

[pks20_4,loc20_4,w20_4,p20_4] =findpeaks(Y20_4,X4)

pks20_4 = 274.4520

loc20_4 = 3.5596

w20_4 = 0.1954

p20_4 = 273.9922

[pks20_4,loc20_4ix,w20_4ix,p20_4] =findpeaks(Y20_4)

pks20_4 = 274.4520

loc20_4ix = 491

w20_4ix = 111.0153

p20_4 = 273.9922

ixv = find(diff(sign(Y20_4-105)))

ixv = 1×2

421 572

PkArea = trapz(X4(ixv(1):ixv(2)), Y20_4(ixv(1):ixv(2)))

PkArea = 50.5773

figure
plot(X4, Y20_4)
grid
hold on
xp = X4(ixv(1):ixv(2));
yp = Y20_4(ixv(1):ixv(2));
patch([xp flip(xp)], [yp zeros(size(yp))], [ 1 1 1]*0.75, 'FaceAlpha',0.5)
hold off
xline(X4(ixv), '--r')
text(mean(X4(ixv(1):ixv(2))), 105, sprintf('Area = %.3f', PkArea), 'Horiz','center')

.

[iqra kiran]

Thank you, can you please elaborate how you have used "105" in this line ixv = find(diff(sign(Y20_4-105))), from where it comes

[Star Strider]

My pleasure!

The ‘105’ value is the approximate y-value where I defined the peak as starting and ending, by looking at it. The ‘ixv’ values are the approximate beginning and ending indices of the peak. The gray region describes the area being calculated (from zero to the peak in that region).

[iqra kiran]

Okay, I have various dataset on which i have to implement this code, this way is'nt it will take time to find the exact range.

how about this way, can you please share you thoughts on this code, results are slight different in these two methods:

[pks20_4,loc20_4,w20_4,p20_4] =findpeaks(Y20_4,X4);

x4=loc20_4(1)-w20_4(1)/2;

x_4=loc20_4(1)+w20_4(1)/2;

area_4 = trapz([x4,x_4], [pks20_4(1), pks20_4(1)]);

[Star Strider]

The default value for the width id ‘half prominence’, so if you want full-width-half-maximum, you would need to specify 'WidthReference','halfheight'.

My analysis —

LD1 = load('X4.mat');
X4 = LD1.X4;
LD2 = load('Y20_4.mat');
Y20_4 = LD2.Y20_4;

[pks20_4,loc20_4,w20_4,p20_4] =findpeaks(Y20_4,X4)

pks20_4 = 274.4520

loc20_4 = 3.5596

w20_4 = 0.1954

p20_4 = 273.9922

[pks20_4,loc20_4ix,w20_4ix,p20_4] =findpeaks(Y20_4)

pks20_4 = 274.4520

loc20_4ix = 491

w20_4ix = 111.0153

p20_4 = 273.9922

ixv = find(diff(sign(Y20_4-105)))

ixv = 1×2

421 572

PkArea = trapz(X4(ixv(1):ixv(2)), Y20_4(ixv(1):ixv(2)))

PkArea = 50.5773

figure
plot(X4, Y20_4)
grid
hold on
xp = X4(ixv(1):ixv(2));
yp = Y20_4(ixv(1):ixv(2));
patch([xp flip(xp)], [yp zeros(size(yp))], [ 1 1 1]*0.75, 'FaceAlpha',0.5)
hold off
xline(X4(ixv), '--r')
text(mean(X4(ixv(1):ixv(2))), 105, sprintf('Area = %.3f', PkArea), 'Horiz','center')

[pks20_4,loc20_4,w20_4,p20_4] =findpeaks(Y20_4,X4);
x4=loc20_4(1)-w20_4(1)/2

x4 = 3.4619

x_4=loc20_4(1)+w20_4(1)/2

x_4 = 3.6572

area_4 = trapz([x4,x_4], [pks20_4(1), pks20_4(1)])

area_4 = 53.6179

ixv = find(diff(sign(Y20_4-p20_4/2)))                               % Find 'Half Prominence' Indices

ixv = 1×2

434 546

ixr = [1 loc20_4ix numel(X4)];
for k = 1:numel(ixr)-1
    idxrng = ixr(k) : ixr(k+1);
    hp(k) = interp1(Y20_4(idxrng), X4(idxrng), p20_4/2);            % Calculate Half-Prominence 'X4' Values
end

hp = 3.4607

hp = 1×2

3.4607 3.6568

HalfPromWidthIdx = X4(ixv(2)) - X4(ixv(1))                          % Half-Prominence Width: Indices

HalfPromWidthIdx = 0.1971

HalfPromWidthItp = hp(2) - hp(1)                                    % Half-Prominence Width: Interpolated Values

HalfPromWidthItp = 0.1961

PkArea = trapz(X4(ixv(1):ixv(2)), Y20_4(ixv(1):ixv(2)))             % Half-Prominence Area: IUndices

PkArea = 42.4737

figure
plot(X4, Y20_4)
grid
hold on
xp = X4(ixv(1):ixv(2));
yp = Y20_4(ixv(1):ixv(2));
patch([xp flip(xp)], [yp zeros(size(yp))], [ 1 1 1]*0.75, 'FaceAlpha',0.5)
hold off
xline(X4(ixv), '--r')
text(mean(X4(ixv(1):ixv(2))), 105, sprintf('Area = %.3f', PkArea), 'Horiz','center')

Using my code with the half-prominence indices results in a different (about 20% smaller) area (42.4737) than ‘area_4’ (53.6179).

Use whatever calculation you wish.

.

[iqra kiran]

Thanks alot.

[Star Strider]

As always, my pleasure!