Here is how I would do it:
% Create Data
x = 0:480;
y = 1E4*(0.5+rand(size(x)));
y(400:450) = 1.3E+5*exp(-(x(400:450)-425).^2/50)+y(400:450);
% Generate Statistics & Identify Peak
ysts = [mean(y) median(y); mean(y)-1.96*std(y)/sqrt(length(y)) mean(y)+1.96*std(y)/sqrt(length(y))];
[ypk, pki] = max(y);
pkiv = find(y >= ysts(2,2)); % Find Peak Indices
Iy = cumtrapz(x, y-ysts(1,1)); % Integrate Entire Record
Iydb = [ones(pkiv(1),1) x(1:pkiv(1))']\Iy(1:pkiv(1))';
Iyd = Iy' - [ones(size(x))' x']*Iydb; % Detrend Using Baseline
Ipk = Iy(max(pkiv)) - Iy(min(pkiv)); % Find Peak Area
% Plot Data & Detrended Integral & Selected Statistics
figure(1)
plot(x, y)
hold on
plot(x, Iyd, '-g')
hold off
legend('Data', 'Detrended Integral', 'Location', 'NW')
text(100, 1E+6, sprintf('Peak Area = %23.15E', Ipk))
The idea is reasonably straightforward:
- Identify the peak as y-values greater than the upper 95% confidence interval for the mean of all the data;
- Integrate the entire record using cumtrapz;
- Do a linear regression on all the integrated data up to the first index of the peak;
- Use that regression to detrend the integrated data;
- Use the identified indices of the peak to calculate the approximate integral of the peak.
This approach may not generalise well if most of the baseline information is not at x-values less than the peak, because it depends on that information to detrend the integral.
The plot produces:
