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
chckxy.m
function [x,y,sizey,endslopes] = mychckxy(x,y) %CHCKXY check and adjust input for SPLINE and PCHIP % [X,Y,SIZEY] = CHCKXY(X,Y) checks the data sites X and corresponding data % values Y, making certain that there are exactly as many sites as values, % that no two data sites are the same, removing any data points that involve % NaNs, reordering the sites if necessary to ensure that X is a strictly % increasing row vector and reordering the data values correspondingly, % and reshaping Y if necessary to make sure that it is a matrix, with Y(:,j) % the data value corresponding to the data site X(j), and with SIZEY the % actual dimensions of the given values. % This call to CHCKXY is suitable for PCHIP. % % [X,Y,SIZEY,ENDSLOPES] = CHCKXY(X,Y) also considers the possibility that % there are two more data values than there are data sites. % If there are, then the first and the last data value are removed from Y % and returned separately as ENDSLOPES. Otherwise, an empty ENDSLOPES is % returned. This call to CHCKXY is suitable for SPLINE. % % See also PCHIP, SPLINE. % Copyright 1984-2011 The MathWorks, Inc. % make sure X is a vector: if length(find(size(x)>1))>1 error(message('MATLAB:chckxy:XNotVector')) end % ensure X is real if any(~isreal(x)) error(message('MATLAB:chckxy:XComplex')) end % deal with NaN's among the sites: nanx = find(isnan(x)); if ~isempty(nanx) x(nanx) = []; warning(message('MATLAB:chckxy:nan')) end n=length(x); if n<2 error(message('MATLAB:chckxy:NotEnoughPts')) end % re-sort, if needed, to ensure strictly increasing site sequence: x=x(:).'; dx = diff(x); if any(dx<0), [x,ind] = sort(x); dx = diff(x); else ind=1:n; end if ~all(dx), error(message('MATLAB:chckxy:RepeatedSites')), end % if Y is ND, reshape it to a matrix by combining all dimensions but the last: sizey = size(y); while length(sizey)>2&&sizey(end)==1, sizey(end) = []; end yn = sizey(end); sizey(end)=[]; yd = prod(sizey); if length(sizey)>1 y = reshape(y,yd,yn); else % if Y happens to be a column matrix, change it to the expected row matrix. if yn==1 yn = yd; y = reshape(y,1,yn); yd = 1; sizey = yd; end end % determine whether not-a-knot or clamped end conditions are to be used: nstart = n+length(nanx); if yn==nstart endslopes = []; elseif nargout==4&&yn==nstart+2 endslopes = y(:,[1 n+2]); y(:,[1 n+2])=[]; if any(isnan(endslopes)) error(message('MATLAB:chckxy:EndslopeNaN')) end if any(isinf(endslopes)) error(message('MATLAB:chckxy:EndslopeInf')) end else error(message('MATLAB:chckxy:NumSitesMismatchValues',nstart, yn)) end % deal with NaN's among the values: if ~isempty(nanx) y(:,nanx) = []; end y=y(:,ind); nany = find(sum(isnan(y),1)); if ~isempty(nany) y(:,nany) = []; x(nany) = []; warning(message('MATLAB:chckxy:IgnoreNaN')) n = length(x); if n<2 error(message('MATLAB:chckxy:NotEnoughPts')) end end
ps:說明下,由於這兩個文件都是matlab的工具文件,所以chckxy.m在調用時,改了名字叫做mychckxy.m,相應的文件名字也需要改。
做一個簡單的測試,做一個調用:
clc; clear all; close all; x = 0:10; y = sin(x); xx = 0:.25:10; yy = spline(x,y,xx) plot(x,y,'o',xx,yy);
運行結果:


到此為止,都是准備工作做,下面開始介紹如何在vs中調用spline函數。
1)在matlab中輸入命令 mbuild -setup , 運行結果如下圖所示, 按照提示選擇編譯器 vs2010.
mbuild -setup

然后鍵入:mex -setup 命令,運行結果如下圖所示,按照提示選擇編譯器 vs2010
mex -setup

然后在matlab命令窗口輸入:
mcc -W cpplib:libspline -T link:lib spline.m
可以得到如下圖這些文件:

依然,其中的"libspline.dll"、"libspline.h"和"libspline.lib"這三個文件是我們所需的。

2)打開vs2010建一個控制台應用程序,可以選擇一個空的控制台應用程序。
創建程序之后把第一步中得到的三個文件copy到工程中。

由於我的電腦是win 7 64bit,matlab是64bit,所以應該選擇x64,而不是win32平台。
a)修改平台參數,為x64
生成 ---> 配置管理器


b)配置包含目錄與庫目錄
項目 ----> 屬性 ----> vc++目錄

包含目錄:
D:\Program Files\MATLAB\R2012a\extern\include


庫目錄:
D:\Program Files\MATLAB\R2012a\extern\lib\win64\microsoft


c)配置附加依賴項 , 這里根據項目的不同,依賴的文件不同,這里測試依賴的是"mclmcrrt.lib"和"libspline.lib"這兩個lib,第一是庫lib,第二個是我們生成的lib.文件。所依賴的lib文件在庫目錄已經說明了,
路徑為:D:\Program Files\MATLAB\R2012a\extern\lib\win64\microsoft下。
這里有兩種解決方案,第一種在vs中配置。
項目 ----> 屬性 ----> 連接器 ----> 輸入


第二種方法是,在文件中直接引入lib文件。

做完以上工作后,我們新建一個主函數作為入口函數,具體測試代碼如下:
#include "libspline.h" //增加頭文件 #include <cmath> #include <iostream> #include <iomanip> using namespace std; #pragma comment(lib,"mclmcrrt.lib") #pragma comment(lib,"libspline.lib") int main() { //初始化lib(必須) if (!libsplineInitialize()) return -1; int i, j; double x[1][11], y[1][11]; for(i=0; i<11; i++) { x[0][i] = i; y[0][i] = sin(x[0][i]); } double xx[1][41]; for(i=0; i<41; i++) xx[0][i] = i*0.25; double yy[1][41]; mwArray mwX(1,11,mxDOUBLE_CLASS); mwArray mwY(1,11,mxDOUBLE_CLASS); mwArray mwXX(1,41,mxDOUBLE_CLASS); mwArray mwYY(1,41,mxDOUBLE_CLASS); mwX.SetData(*x, 11); mwY.SetData(*y, 11); mwXX.SetData(*xx, 41); mwYY.SetData(*yy, 41); spline(1, mwYY, mwX, mwY, mwXX); //調用spline cout<<"yy = "<<endl; i = 0; for(j = 0; j < 41; j++) { //Get第一個參數表示用1個下標訪問元素,j+1是列號(MATLAB下標從1開始,而C++從0開始,故做+1操作) yy[0][j] = mwYY.Get(1,j+1); cout<<setprecision(4)<<right<<setw(10)<<yy[0][j]; i++; if(i%7 == 0) cout<<endl; //換行 } cout<<endl; //終止調用 libsplineTerminate(); return 0; }
運行結果如圖:

比較這個結果與最開始我們測試matlab運行的結果,測試通過。matlab配置完成。
ps說明:配置過程中遇到的問題:
配置時經常遇到 LINK2019的錯誤。這種錯誤就是典型的lib缺失導入的問題。
main.obj : error LNK2019: 無法解析的外部符號 mclGetMatrix_proxy,該符號在函數 "public: __cdecl mwArray::mwArray(unsigned __int64,unsigned __int64,enum mxClassID,enum mxComplexity)" (??0mwArray@@QEAA@_K0W4mxClassID@@W4mxComplexity@@@Z) 中被引用 1>main.obj : error LNK2019: 無法解析的外部符號 mclcppGetLastError_proxy,該符號在函數 "public: static void __cdecl mwException::raise_error(void)" (?raise_error@mwException@@SAXXZ) 中被引用 1>main.obj : error LNK2019: 無法解析的外部符號 mclcppCreateError_proxy,該符號在函數 "public: __cdecl mwException::mwException(void)" (??0mwException@@QEAA@XZ) 中被引用 1>main.obj : error LNK2019: 無法解析的外部符號 ref_count_obj_addref_proxy,該符號在函數 "public: __cdecl mwException::mwException(class mwException const &)" (??0mwException@@QEAA@AEBV0@@Z) 中被引用 1>main.obj : error LNK2019: 無法解析的外部符號 ref_count_obj_release_proxy,該符號在函數 "public: virtual __cdecl mwException::~mwException(void)" (??1mwException@@UEAA@XZ) 中被引用 1>main.obj : error LNK2019: 無法解析的外部符號 error_info_get_message_proxy,該符號在函數 "public: virtual char const * __cdecl mwException::what(void)const " (?what@mwException@@UEBAPEBDXZ) 中被引用 1>main.obj : error LNK2019: 無法解析的外部符號 array_ref_getV_int_proxy,該符號在函數 "public: class mwArray __cdecl mwArray::GetPromoted(unsigned __int64,...)" (?GetPromoted@mwArray@@QEAA?AV1@_KZZ) 中被引用 1>main.obj : error LNK2019: 無法解析的外部符號 array_ref_set_numeric_mxDouble_proxy,該符號在函數 "public: void __cdecl mwArray::SetData(double *,unsigned __int64)" (?SetData@mwArray@@QEAAXPEAN_K@Z) 中被引用 1>main.obj : error LNK2019: 無法解析的外部符號 array_ref_get_numeric_mxDouble_proxy,該符號在函數 "public: __cdecl mwArray::operator double(void)const " (??BmwArray@@QEBANXZ) 中被引用
這里是因為缺少:mclmcrrt.lib
#pragma comment(lib,"mclmcrrt.lib")
即可解決。
