1、背景
众所周知,matlab在处理矩阵、数学计算、计算机仿真、图像处理等方面有着 c c++无可比拟的优势,但是做成系统供使用时,又显得过于粗糙,为了使用起来高大上,计算起来有简单,方便。无疑,c++ 与matlab混合编程将会使非常靠谱的选择。
这里暂且不论所谓的matlab效率低,c/c++效率高的问题,自我感觉,以我目前编码的功底,所编写的代码的效率远远不及matlab提供的代码的效率。除非你是大牛,或者你是人云亦云,所以能用matlab混合c++编码还是很不错的选择,话不多说,我们开始讨论正题。
2、我使用的版本是matlab2012与vs2010混合编程的。
软件的下载这里就不多说了,我相信看这篇教程的你,这两个软件已经安装的妥妥当当的了。
这里我选用网上常用来做例子的matlab代码做测试,spline.m,该文件位于
D:\Program Files\MATLAB\R2012a\toolbox\matlab\polyfun
当然该文件中依赖调用另一个文件chckxy.m,该文件也在这条路径下。找到后复制到matlab的工作目录下。
这里为了方便提供两个文件的代码:
spline.m
function output = spline(x,y,xx) %SPLINE Cubic spline data interpolation. % PP = SPLINE(X,Y) provides the piecewise polynomial form of the % cubic spline interpolant to the data values Y at the data sites X, % for use with the evaluator PPVAL and the spline utility UNMKPP. % X must be a vector. % If Y is a vector, then Y(j) is taken as the value to be matched at X(j), % hence Y must be of the same length as X -- see below for an exception % to this. % If Y is a matrix or ND array, then Y(:,...,:,j) is taken as the value to % be matched at X(j), hence the last dimension of Y must equal length(X) -- % see below for an exception to this. % % YY = SPLINE(X,Y,XX) is the same as YY = PPVAL(SPLINE(X,Y),XX), thus % providing, in YY, the values of the interpolant at XX. For information % regarding the size of YY see PPVAL. % % Ordinarily, the not-a-knot end conditions are used. However, if Y contains % two more values than X has entries, then the first and last value in Y are % used as the endslopes for the cubic spline. If Y is a vector, this % means: % f(X) = Y(2:end-1), Df(min(X))=Y(1), Df(max(X))=Y(end). % If Y is a matrix or N-D array with SIZE(Y,N) equal to LENGTH(X)+2, then % f(X(j)) matches the value Y(:,...,:,j+1) for j=1:LENGTH(X), then % Df(min(X)) matches Y(:,:,...:,1) and Df(max(X)) matches Y(:,:,...:,end). % % Example: % This generates a sine-like spline curve and samples it over a finer mesh: % x = 0:10; y = sin(x); % xx = 0:.25:10; % yy = spline(x,y,xx); % plot(x,y,'o',xx,yy) % % Example: % This illustrates the use of clamped or complete spline interpolation where % end slopes are prescribed. In this example, zero slopes at the ends of an % interpolant to the values of a certain distribution are enforced: % x = -4:4; y = [0 .15 1.12 2.36 2.36 1.46 .49 .06 0]; % cs = spline(x,[0 y 0]); % xx = linspace(-4,4,101); % plot(x,y,'o',xx,ppval(cs,xx),'-'); % % Class support for inputs x, y, xx: % float: double, single % % See also INTERP1, PCHIP, PPVAL, MKPP, UNMKPP. % Carl de Boor 7-2-86 % Copyright 1984-2010 The MathWorks, Inc. % $Revision: 5.18.4.6 $ $Date: 2010/09/02 13:36:29 $ % Check that data are acceptable and, if not, try to adjust them appropriately [x,y,sizey,endslopes] = mychckxy(x,y); n = length(x); yd = prod(sizey); % Generate the cubic spline interpolant in ppform dd = ones(yd,1); dx = diff(x); divdif = diff(y,[],2)./dx(dd,:); if n==2 if isempty(endslopes) % the interpolant is a straight line pp=mkpp(x,[divdif y(:,1)],sizey); else % the interpolant is the cubic Hermite polynomial pp = pwch(x,y,endslopes,dx,divdif); pp.dim = sizey; end elseif n==3&&isempty(endslopes) % the interpolant is a parabola y(:,2:3)=divdif; y(:,3)=diff(divdif')'/(x(3)-x(1)); y(:,2)=y(:,2)-y(:,3)*dx(1); pp = mkpp(x([1,3]),y(:,[3 2 1]),sizey); else % set up the sparse, tridiagonal, linear system b = ?*c for the slopes b=zeros(yd,n); b(:,2:n-1)=3*(dx(dd,2:n-1).*divdif(:,1:n-2)+dx(dd,1:n-2).*divdif(:,2:n-1)); if isempty(endslopes) x31=x(3)-x(1);xn=x(n)-x(n-2); b(:,1)=((dx(1)+2*x31)*dx(2)*divdif(:,1)+dx(1)^2*divdif(:,2))/x31; b(:,n)=... (dx(n-1)^2*divdif(:,n-2)+(2*xn+dx(n-1))*dx(n-2)*divdif(:,n-1))/xn; else x31 = 0; xn = 0; b(:,[1 n]) = dx(dd,[2 n-2]).*endslopes; end dxt = dx(:); c = spdiags([ [x31;dxt(1:n-2);0] ... [dxt(2);2*(dxt(2:n-1)+dxt(1:n-2));dxt(n-2)] ... [0;dxt(2:n-1);xn] ],[-1 0 1],n,n); % sparse linear equation solution for the slopes mmdflag = spparms('autommd'); spparms('autommd',0); s=b/c; spparms('autommd',mmdflag); % construct piecewise cubic Hermite interpolant % to values and computed slopes pp = pwch(x,y,s,dx,divdif); pp.dim = sizey; end if nargin==2, output = pp; else output = ppval(pp,xx); end