Don't like ads? PRO users don't see any ads ;-)
Guest

Untitled

By: a guest on May 28th, 2012  |  syntax: None  |  size: 2.35 KB  |  hits: 14  |  expires: Never
download  |  raw  |  embed  |  report abuse  |  print
Text below is selected. Please press Ctrl+C to copy to your clipboard. (⌘+C on Mac)
  1. function [xx, yy] = mySplinePlot(x, y, accuracyFactor)
  2.  
  3. % Sort the points by increasing x-coordinates
  4. [x ind] = sort(x);
  5. y = y(ind);
  6.  
  7. % Plot the original points
  8. subplot(2,2,1)
  9. plot(x,y,'ro')
  10. hold on
  11. title('mySplinePlot')
  12.  
  13.  
  14.  
  15. % Find the number of points based on the max difference between x-values
  16. % and create an xx vector of this number of points, evenly spaced. The max
  17. % difference between adjacent points should be rounded up, as should the
  18. % result after applying the accuracy factor.
  19. pointsPerInterval = round((round(max(diff(x)))) .* accuracyFactor);
  20. totalPoints = (length(x)-1)*pointsPerInterval + length(x);
  21. xx = linspace(x(1), x(end), totalPoints);
  22.  
  23.  
  24. % Interpolate to find a yy vector using the Lagrange Interpolation
  25. % Polynomial. Initialize yy at the right size for improved efficiency.
  26. yy = zeros(1,totalPoints);
  27. for i = 1:totalPoints
  28. yy(i) = lagrangeInterp(x, y, xx(i));
  29. end
  30.  
  31. % Plot the vectors just calculated.
  32. plot(xx, yy)
  33. hold off
  34.  
  35.  
  36.  
  37. % Call Matlab's spline and plot it as described in the problem statement
  38. splineYY = spline(x, y, xx);
  39. subplot(2,2,2)
  40. hold on
  41. plot(x,y,'ro')
  42. plot(xx, splineYY)
  43. title('spline')
  44. hold off
  45.  
  46. % Make a plot as described in the problem statement and as shown in the
  47. % solution figures comparing mySplinePlot with spline
  48. subplot(2,1,2)
  49. plot(xx, yy, 'r')
  50. hold on
  51. plot(xx, splineYY, 'k--')
  52. legend('mySplinePlot','spline')
  53. title('mySplinePlot vs. spline')
  54. plot(x,y,'ro')
  55. hold off
  56.  
  57. end
  58.  
  59.  
  60. % This function performs the actual Lagrange Interpolation, returning an
  61. % interpolated y-value for a single given x-value. Do NOT change this
  62. % function header.
  63. function interpY = lagrangeInterp(x, y, interpX)
  64.  
  65. % Initialize interpY and determine the number of terms in the product The
  66. % new interpolated Y value will then be increased based on the requirements
  67. % of the Lagrange Interpolation.
  68. interpY = 0;
  69. n = length(x);
  70.  
  71. % Make an estimation for each known x-coordinate using the formula given.
  72. for i = 1:n
  73.  
  74. % If this makes more sense to you in a for loop, you can include a
  75. % nested for loop here. This is where the actual product from the
  76. % Lagrange Interpolation formula comes in.
  77. product = prod(interpX - [x(1:i-1), x(i+1:end)])/prod(x(i) - [x(1:i-1), x(i+1:end)]);
  78.  
  79. % Update the total for interpY based on the value found for the ith
  80. % coordinate
  81. interpY = interpY + product * y(i);
  82. end
  83. end