Assignment 2 Scripts
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% Problem A %%%%%%%%%%%%%%%%%%%%%% % Script ShowBigDipper clear % Big Dipper Evolution % The Big Dipper is a constellation consisting of % seven stars which we index as follows: % % % 2 * 3* 7 * % 1 * 4 * % % 6 * % 5 * % The location of these stars varies with time. We have data % for times t(1),...,t(5) where t = [-50; -25; 0; 25; 50]; % Think of t as (-50000BC, -25000BC, Now, 25000AD, 50000AD) % The (x,y) coordinate of star j at time t(i) is specified by % (x(i,j),y(i,j)) where x =[ 15.5 29.0 36.2 45.2 49.5 62.0 55.2;... 16.0 28.2 34.9 44.6 49.2 61.5 56.5;... 17.8 27.9 34.3 44.1 48.6 61.0 58.9;... 18.5 27.5 33.9 43.8 48.5 60.5 61.2;... 20.3 26.2 33.2 42.9 47.8 59.7 62.5]; y = [ 0.0 5.0 5.1 5.5 0.0 6.7 14.8 ;... -0.4 5.3 4.3 5.2 0.0 4.6 14.3 ;... 0.5 5.5 4.3 5.2 0.0 4.1 13.1 ;... 0.4 5.4 4.3 5.7 0.0 4.4 12.7 ;... -1.2 5.5 4.4 6.2 0.0 5.0 13.2]; % Here is a 5-frame "movie" of the evolution... bigscreen % An NCM function that makes the figure window big starsize = 18; % Might want to play with this value. fps = 1; % frames per second close all for k=1:5 hold on fill([0 80 80 0 0],[-30 -30 50 50 -30],'k') plot(x(k,:),y(k,:),'r') plot(x(k,:),y(k,:),'*y','markersize',starsize) text(30,40,sprintf('Time = %6.0f',1000*t(k)),'color','y','fontsize',24) hold off axis([0 80 -30 50]) axis equal off pause(2) % Take a moment to look at the frame M(k) = getframe; end close all movie(M,3,1) % Play the movie % A2A % Shows Big Dipper Evolution movie bigscreen fps = 1; % frames per second r = 3; % movie shown this number of times M = BigDipperMovie(linspace(-50,100,16)'); shg movie(M,r,fps) %%%%%%%%%%%%%%%%%%%%% Problem B %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% function ShowMultipleInterpN(n,m) % Illustrates MultipleInterpN % n and m are integers that satisfy 1 < = m < = n % % m-point polynomial interpolants of sin(x) across % [0,2pi] are displayed. n shouldn't be too much bigger % than 20. x = linspace(0,2*pi,n)'; y = sin(x); C = MultipleInterpN(x,y,m); x0 = linspace(0,2*pi)'; y0 = sin(x0); [p,q] = size(C); for k=1:q xsub = x(k:k+m-1); pVals = HornerN(C(:,k),xsub,x0); subplot(2,1,1) plot(x0,y0,x0,pVals,xsub,sin(xsub),'o',x,sin(x),'*k') title(['Interpolation point indices = ' sprintf('%1d ',k:k+m-1)]) axis([0 2*pi -3 3]) subplot(2,1,2) semilogy(x0,abs(pVals-y0)+eps) axis([0 2*pi 10^-6 10^1]) title('Approximate Absolute Error') shg pause end %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% problem C %%%%%%%%%%%%%%%%%%%% % Script A2C % Check for correctness.. A = linspace(0,360,13)'; D = [ 408 89 -66 10 338 807 1238 1511 1583 1462 1183 804 408]'; Avals = linspace(0,360,200)'; F = newCSInterp(D(1:12)); Fvals = CSeval(F,360,Avals); plot(Avals,Fvals,A,D,'o') axis([-10 370 -200 1700]) set(gca,'xTick',linspace(0,360,13)) xlabel('Ascension (Degrees)') ylabel('Declination (minutes)') % Check efficiency nRepeat = 1; % Number of instances that are timed % Make big enough so clock "granularity" is not an issue... clc disp(' n CSInterp Time/ newCSInterp time') disp('----------------------------------------') for n = 100:50:500 f = randn(n,1); % Benchmark the original implementation on page SCMV p. 101 tic for k=1:nRepeat F = CSInterp(f); end t1 = toc; tic for k=1:nRepeat F = newCSinterp(f); end t2 = toc; disp(sprintf(' %3d %6.3f',n,t1/t2)) end