- Woops, in case anyone is using my attached files: the CSV-file that the .m file is loading is not the one I put up. It shouldn't matter though.
Finding area under peaks for drifting signal
6 次查看(过去 30 天)
显示 更早的评论
EDITED, see bottom, (regarding numerical integration)
Dear all,
How can I find the area under the "larger" peaks, as seen below?
![](https://www.mathworks.com/matlabcentral/answers/uploaded_files/169345/image.jpeg)
My idea is to:
1) locate where the peaks are (I have tried using the findpeaks() function)
2) add or subtract the minimum value of the particular peak to get a baseline for each peak at y = 0.
3) use some sort of trapezoidal integration (I saw that MATLAB has a built-in function that does this).
When using the findpeaks() function, I'm having a really hard time locating "interesting" peaks. For instance, the first "big" peak has some noise in it and then findpeaks() give me a lot of peaks on it - and then there are other small peaks I'm not interested in.
I hope someone can help!
(If you need it, my data is found in the attached CSV file - I'm only using the first two columns, "Time" and "Channel1V". I have 100 CSV files like the uploaded one, so I really need a smart code ;) )
EDIT:
Q: How to find area under peak using trapz()?
So, using the arguments 'MinPeakDistance' and 'MinPeakProminence' I have been able to figure out which peaks to use. I still need to find the area under the graphs; so far, I've been using the very rough estimate Prominence times Width (where width is the width of the peak taken at the half of the prominence). It should for most purposes be alright, but I'd really appreaciate if anyone has got a good idea of how to numerically integrate this.
11 个评论
Star Strider
2018-1-25
My pleasure.
I can help with the signal processing. You posted one of your signals, so I will download it, and experiment with it.
采纳的回答
Star Strider
2018-1-25
I would use the Signal Processing Toolbox medfilt1() function to estimnate the base line drift, that you can then remove. (None of my other approaches, such as frequency-selective filtering or the Savitzky-Golay filter, gave satisfactory results.)
Try this to remove the baseline drift:
[D,S,R] = xlsread('420B1.1.csv');
T = D(:,1);
Ch1 = D(:,2);
Ch2 = D(:,3);
Ts = mean(diff(T));
L = length(T);
MF = medfilt1(Ch1, 750); % Median Filter
Ch1F = Ch1-MF; % Correct Baseline
figure
plot(T,Ch1F)
grid
I leave the peak detection to you. Here is one possibility:
dCh1F = gradient(Ch1F, Ts); % Derivative
[pkp,locp] = findpeaks(dCh1F, 'MinPeakHeight',0.04, 'MinPeakDistance',500); % Positive Derivative Peaks, ‘locp’ Are Indices
[pkn,locn] = findpeaks(-dCh1F, 'MinPeakHeight',0.03, 'MinPeakDistance',500); % Negative Derivative Peaks, ‘locn’ Are Indices
figure
plot(T,Ch1F)
hold on
plot(T([locp; locn]),Ch1F([locp; locn]), '+r')
hold off
grid
You will have to experiment to get the result you want. This may not be robust to all your data, so you may have to adjust it for each file.
2 个评论
更多回答(1 个)
AJ Geiger
2018-1-24
1. instead of using the function maxpeak use
[~,I] = max(f);
to obtain the index of the max peak
2. Nice! Now you know where the max peak is. The next step is to set up your boundaries, so you know where to start and stop your integration. There are several algorithms you can implement. the easiest would be simulated annealing. The algorithm goes as follows.
% Note: E is your data vector
% for the alg I ~= 1
if I == 1 % check condition
I = 2;
end
%
%
% Preform Alg
for t = 0:length(E)-1
% check if next point is
idx = I+t;
Delta_E = E(I+t) - E(I-1+t);
if E > 0, continue, end
%
% Note T: is your cooling function which
% Converges to 0 as t increases. now this
% is basic but T could be as simple as
% T = 1/t or T = 1/(t^2)
% now for tuning part to be complete you
% should say:
% if T == 0, break, end
% however only use it if you develop T
% properly.
%
T = 1/(t^2); % Temperature Function
if exp(Delta_E/T) < Thresh, break, end
end
Now you have found the Right-hand side boundary. For now, you can
3. To be complete, idx_2 = idx;
4. use the function E_New = flip(E);
5. solve for I_New = I - length(E) + 1;
6. then run alg two on E_new and I_New
7. To be complete, idx_1 = idx;
8. get index vector II you wish to integrate over
II = idx_1:idx_2
9. Then solve for the area using II
Area = trapz(E(II));
10. Grab a beer if of legal age and does not effect believes (Optional)
11. Finish beer
12. Explore other integration methods.
另请参阅
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!