RoundTrip - function RoundTrip Program to add noise to a...

Info iconThis preview shows pages 1–3. Sign up to view the full content.

View Full Document Right Arrow Icon
Sheet1 Page 1 function RoundTrip % Program to add noise to a sinusoidal test function and take two % derivatives, and then use numerical integration method to take % two integrals to see how well we recover the original function. % Define sinusoidal test function npts = 21 t = linspace(0,1,npts)' a = 1 b = 2 tf = 1 w = 2*pi/tf q = a*sin(w*t)+b % Add noise to original data noisemag = 0.05 qnoisy = q+noisemag*(2*rand(npts,1)-1) plot(t,q,'k-','LineWidth',2) hold on plot(t,qnoisy,'ko','LineWidth',2) set(gca,'FontSize',14) legend('Original','Noisy') pause % Calculate first derivative of original data analytically and numerically % and of noisy data numerically. Use forward difference for first % derivative calculations. qp = w*a*cos(w*t) h = t(2,1)-t(1,1) qporig = zeros(npts,1) qpnoisy = zeros(npts,1) for i = 1:npts if i == npts % Use second data point as first sample point by assuming % periodicity qporig(i,1) = (q(2,1)-q(i,1))/h qpnoisy(i,1) = (qnoisy(2,1)-qnoisy(i,1))/h else qporig(i,1) = (q(i+1,1)-q(i,1))/h qpnoisy(i,1) = (qnoisy(i+1,1)-qnoisy(i,1))/h end end hold off plot(t,qp,'k-','LineWidth',2) hold on plot(t,qporig,'ko','LineWidth',2) plot(t,qpnoisy,'kx','LineWidth',2) set(gca,'FontSize',14)
Background image of page 1

Info iconThis preview has intentionally blurred sections. Sign up to view the full version.

View Full DocumentRight Arrow Icon
Sheet1 Page 2 legend('Analytical','Numerical orig','Numerical noisy') title('Forward difference - first derivative') % Note that we have two sources of error: % 1) Discretization with finite step size h (and note the systematic nature % of the error), and % 2) Noise in the data. pause % Now repeat the first derivative calculations using central differencing for i = 1:npts if i == 1 % Use next to last data point as second sample point by assuming % periodicity qporig(i,1) = (q(i+1,1)-q(npts-1,1))/(2*h) qpnoisy(i,1) = (qnoisy(i+1,1)-qnoisy(npts-1,1))/(2*h)
Background image of page 2
Image of page 3
This is the end of the preview. Sign up to access the rest of the document.

This note was uploaded on 09/05/2011 for the course EGM 3344 taught by Professor Raphaelhaftka during the Spring '09 term at University of Florida.

Page1 / 5

RoundTrip - function RoundTrip Program to add noise to a...

This preview shows document pages 1 - 3. Sign up to view the full document.

View Full Document Right Arrow Icon
Ask a homework question - tutors are online