필터 지우기
필터 지우기

How to interpolate midpoints in a curve

조회 수: 5 (최근 30일)
Proman
Proman 2020년 9월 8일
댓글: Proman 2020년 9월 9일
Hello to everyone and greetings
I have a curve that I intend to find its midpoints (points between each succesive minimum and maximum) and connect them to each other to plot an envelope-wise curve like the following (the green dotted-dashed curve is the envelope I try to realize):
the code I used to plot this is the trivial approach (to my opinion):
clc
clear
close all
load Data1
figure
plot(t,R66)
xlabel('time (ps)')
ylabel('\rho_{66}')
[pks,max_loc] = findpeaks(R66,t);
invert_R66 = max(R66(:)) - R66;
[pksm,min_loc] = findpeaks(invert_R66,t);
pksm = max(R66(:)) - pksm;
xmid = (max_loc(2 : end) + min_loc) ./ 2;
xmid_p = (min(R66) + max_loc(1)) / 2;
ymid = interp1(t,R66,xmid);
ymid_p = interp1(t,R66,xmid_p);
xmid = [xmid_p xmid'];
ymid = [ymid_p ymid'];
hold on
plot(xmid,ymid,'r')
axis([0 700 0 0.045])
legend('Population','Restoration')
It quite worked well for every curve I gave as an input but the following curve (which its data is attached as Data1.mat) is not responding well as you can see
I somehow believe that there must be a problem with interpolation but after all, I can not figure out where I am going wrong. It would be kind of you to help me at this.
Thanks in advance for your time given on my question.
  댓글 수: 2
Matt J
Matt J 2020년 9월 8일
편집: Matt J 2020년 9월 8일
The Data1.mat file contains several data sets. Which is the one that doesn't work?
Proman
Proman 2020년 9월 9일
R66 is the one that goes wrong. Already indicated in the code

댓글을 달려면 로그인하십시오.

채택된 답변

Matt J
Matt J 2020년 9월 9일
편집: Matt J 2020년 9월 9일
The code is working fine. You think there is only one peak in the interval 200<=t<=300, but if you zoom in, you will see that there are three of them squished very close together. The same for the other intervals.
  댓글 수: 4
Matt J
Matt J 2020년 9월 9일
편집: Matt J 2020년 9월 9일
A Savitzky-Golay pre-filter seems to help a lot.
load('Data1.mat')
close all
figure
procfun(R66,t)
function procfun(R,t)
w=101; %window width
Rsg = sgolayfilt(R,2,w);
Rsg(1:w)=R(1:w);
D=diff(sign(diff([Rsg(1);Rsg])))~=0;
tp=t(D);
tq=(tp(1:end-1)+tp(2:end))/2;
Rq=interp1(t,R,tq);
Rp=interp1(t,R,tp);
pp=interp1(tq,Rq,'linear','pp');
hold on
plot(t,R,tp,Rp,'rx');
fplot(@(t)ppval(pp,t),[tq(1),t(end)])
xlabel('time (ps)')
ylabel('\rho_{66}')
hold off
end
Proman
Proman 2020년 9월 9일
What a miracle this filter is *_*
Thanks again for your care and extremely productive help so far my friend!

댓글을 달려면 로그인하십시오.

추가 답변 (0개)

카테고리

Help CenterFile Exchange에서 Spectral Measurements에 대해 자세히 알아보기

제품


릴리스

R2019b

Community Treasure Hunt

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

Start Hunting!

Translated by