Smoothing piecewise sin functions
16 views (last 30 days)
Show older comments
I wrote a function that pieces together two different sine functions. The positive curves are always sin(x) but the negative curves are sin(mult.*x).
I tried spaps/csaps/etc and spaps with tol=1 seems to give me the smoothest result, except the amplitude of the curve decreases from the desired y=-1:1. I tried adjusting the tol to different values but that didn't quite work. Is there some other function/parameter I could use to preserve both the smoothness at tol=1 and the amplitude?
The code is below:
function [t, values]=graphic(mult);
clf
x=0:pi./32:4.*pi;
y=[sin(x).' x.'];
%y=y(y>=0);
for m=1:length(x)
if y(m) >=0
y(m)=y(m);
else
y(m)=0;
end
end
xx=0:pi./32:4.*pi;
yy=[sin(mult.*xx).' xx.'];
%yy=yy(yy<0);
for n=1:length(xx)
if yy(n) < 0
yy(n)=yy(n);
else
yy(n)=0;
end
end
Y=([yy y]);
SET1=find(Y(:,3)~=0);
SET2=find(Y(:,1)~=0);
tol=1e-15;
ZERO1=find(Y(SET1,3)<tol);
ZERO2=find(Y(SET2(1):end,1)==0,1);
Inter=[Y(SET1(1:ZERO1(1)),3); Y(SET2(1:ZERO2(1)),1)];
INTER=[0; Inter; Inter];
t=1:length(INTER);
%plot(t,INTER);
%y2=smooth(INTER,'sgolay')
% p=polyfit(t,INTER.',6);
% y2=polyval(p,t);
%plot(t,y2);
[sp,values]=spaps(t,INTER,1);
fnplt(sp);
set(gca,'XTick',0:25:length(INTER));
end
2 Comments
dpb
on 23 Aug 2013
Edited: dpb
on 24 Aug 2013
Don't think your above code works or if does run doesn't do what you think/want...
"The Matlab way"... :)
x=0:pi/32:4*pi;
y=sin(x)';
y(y<0)=0; % above you've already built a 2-col matrix but
y=[y x']; % address it only w/ one subscript.
Not sure what the purpose of the x in the y array is for???
Ditto comments above for next section plus your 'xx' is identical to 'x'; why not just use the copy you already have?
yy=sin(mult.*x)'; % do the fixup before the concatenation
yy(yy>=0)=0;
yy=[yy x']; % again, why the concatenation into yy
why don't you just write
x=[0:pi/32:4*pi]'; % Start w/ a column if that's what want...
y=sin(x); % then y is already column, too.
iy=y<0; % logical vector of the y<0 locations
y(iy)=0; % clear them
y(~iy)=sin(mult*x(~iy)); % fill in the alternate locations
What's mult? Not clear how to answer the rest of the question w/o that minor detail...
Answers (2)
Image Analyst
on 23 Aug 2013
Your code doesn't run, and you didn't upload a diagram, so it's hard to imagine what it looks like. But I imagine stretches of sine waves with different periods (frequencies) and where they meet you are going to have a huge step in intensity.
To smooth the steps out, you can use Fourier filtering (or use conv2()), or you can use a Savitzky-Golay filter. A sine can be pretty well fit by a 3rd or 5th order polynomial in the short run, so if you have the Signal Processing Toolbox, you can use sgolay() with an order 3 or 5 and it should do a pretty good job of leaving the sine wave parts along but smoothing over the abrupt steps where the sines of different frequencies meet.
9 Comments
dpb
on 24 Aug 2013
Edited: dpb
on 24 Aug 2013
Instead of trying to mung on the result of mismatching frequencies caused by the discontinuous mult value as you have, redefine it to be a continuous function w/ a rapid change.
A sigmoid or the like would work to remove the discontinuity but still allow a rapid change from one to the other frequency and vice versa.
ADDENDUM: Example for first cycle of switch from wx=1x to wx=2x as a function of changing the intensity of the sigmoid. Try the following from command line...
x=[0:pi/128:2*pi]'; % need better resolution than your original...
N=[0.1 0.5 1 2 5 10]'; % a set of weights for the sigmoid--read code
M=2; % the target multiplier value
s=zeros(length(x),length(N)); % preallocate so can plot all together later
for i=1:length(N)
n=N(i);
w=(M-1)*(tanh(n*(x-pi))+1)/2+1; % the omega weighted vector
s(:,i)=sin(w.*x); % and the sin resulting
end
plot(x,s)
legend(num2str(N))
You'll observe that as you make the transition sharper and sharper the resulting discontinuity approaches that of the limiting case--I had envisioned you could get by "more better" than it actually turns out, but certainly the way to generate a smooth transition is something on this order rather than trying to fix up the abrupt transition as currently trying to do.
Doing the obverse transition back at 2pi is left as "exercise for the student" -- enjoy :)
OBTW, as did in previous comment, it's worthwhile to plot the weight function and/or just the sigmoid as function of the various parameters to compare the effects w/ the result.
ADDENDUM 2:
On the comment of modifying the pi:2pi regions only, change the offset bias from the above pi to something >pi: that will be the breakpoint. Plotting the sigmoid for various choices of parameters is best way to see what's going on and what affects what...
8 Comments
dpb
on 12 Sep 2013
Couple of ways...
a)
w2=(M-1)*(tanh(n*(x(floor(end/2):end)-3*pi))+1)/2+2;
where now
x=0:4*pi
interval so w2 is covering 2pi-4pi and w is as before except for 0-2pi or x(1:round(end/2))
b) x as above over the extended range, then create w as
w=[w;flipud(w)];
May need to fixup length(w) by an "off-by-one" error if length(x) is odd the length of doubled w will be one too great so use (end-1:-1:1) for the added portion to match the length and not duplicate the middle point.
See Also
Categories
Find more on Spline Postprocessing in Help Center and File Exchange
Products
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!