sba一个通用的稀疏光束法平差的软件包

sba一个通用的稀疏光束法平差的软件包 如果你来到这个页面来寻找一个通用的Levenberg-Marquardt算法的C/C++实现，请看levmar 引言: 本页面是关于sba，一个通用的稀疏光束法平差的C/C++软件包。它基于GNU通用公共许可证GPL分发的。光束法平差(BA)是作为每个基于特征的多视重建视觉算法的最后一步，用来获得最佳的三维结构和运动(如相机矩阵)参数估计。提供初始估计，BA同时精化运动和结构参数，通过最小化观测和预测的图像点之间的投影误差。最小化一般通过Levenberg-Marquardt (LM)算法来辅助完成。然而，由于许多未知的因素作用于最小投影误差，一个通用的LM算法的实现(如MINPACK的 lmder)当应用于BA背景下的定义的最小化问题时，会带来极高的计算代价。 幸运的是，在基本的法方程中不同的三维点和相机参数相互之间影响较小，呈现一种稀疏的块结构(如图)。Sba利用这种稀疏的特性，使用LM算法的简化的稀疏变量来降低计算的复杂度。Sba是通用的，因为它保证了用户对于相机和三维结构的描述参数的定义的完全控制。因此，它事实上可以支持任何多视重建问题的显示和参数化。比如任意投影相机，部分的或完全标定的相机，由固定的三维点进行外方位元素(即姿态)的估计，精化本征参数，等等。用户要想在这类问题中使用sba，只需要提供合适的程序对这些问题和参数来计算估计的图像投影和他们的函数行列式(Jacobian)。用来计算解析的函数行列式可以是手头的代码，或者使用支 持符号微分的工具(如maple)生成的代码，或者通过自动微分技术获得的代码。也可以使用近似的函数行列式，辅之以有限差分的 方法 快递客服问题件处理详细方法山木方法pdf计算方法pdf华与华方法下载八字理论方法下载 。另外，sba包含了检查用户提供的函数行列式的一致性的程序。就我们的知识之所及，sba是第一个并且也是当前独一无二的的软件包，因为他能够不受版权限制以源代码形式放置在任何工程中。 作为sba的效率的一个指标，我们在这里说明，sba的单次测试已经涉及54台相机和5207三维点，产生了24609个图像投影。相应的最小化问题依赖于15999个变量，sba使用非最优的BLAS在 Intel P4@1.8 GHz running Linux机器上大约7秒钟内解决。如果没有BA的稀疏实现，那么这种规模的问题会变得非常棘手。 又一个光束法平差库，由德国斯图加特大学发布。名字很怪，不知道全称是什么。 引言 程序DGAP 实现了光束法平差的摄影测量方法，由 Helmut Schmid and Duane Brown 发明。它 基于图像和目标的几何关系的中心投影，使用最小二乘法。 特点 Camera-/self-/simultan 标定，连同作者 Brown, Ebner and Gruen 建议附加的参数。 两者可选的图像模型:直接线性变换(DLT)和仿射变换。 不同测量位置和/或姿态数据(GPS 支持的空三、直接地理参考)的集成。 精确的计算内外几何参数 测试附加参数的意义 计算协方差(新～)和相关性 新:对分格摄影机图像(frame camera imagery)的扩展的摄影测量模 型。 新:线扫描仪(line scanner)图像的直接地理参考。 新:空三(aerial triangulation，AT)的样例，GPS 支持的空三并且直 接地理参考。 新:ADS-40 线扫描仪图像的直接地理参考样例。 版权/许可 Copyright (C) 2005 Dirk Stallmann This program is free software; you can redistribute it and/or modify it under the terms of the GNU General Public License as published by the Free Software Foundation; either version 2 of the License, or (at your option) any later version. This program is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License for more details. You should have received a copy of the GNU General Public License along with this program; if not, write to the Free Software Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA 02111-1307, USA. 下载地址: tml 依赖 库genlib2 提供了各种子程序以及模板类。 可以有选择的使用LAPACK 和BLAS。LAPACK 是一个Fortran77 库，用来解数字线性代数中 最常见的问题。它又依赖于基本的线性代数子程序BLAS。BLAS 也是一个Fortran77 库，提供 优化的向量和矩阵操作。两个库都是自由软件，并且拥有版权。 LAPACK 和 BLAS 可在netlib 上获得: 为了在MS Windows 2000/XP 下编译DGAP，可以使用Cygnus 的GNU-win32 开发工具包 “Cygwin” ()，版本1.3.6 或者更高。 LAPACK/BLAS 库在Cygwin 下 依然可用。 安装 在UNIX 或者Cygwin 下编译 1、在同一目录下解压压缩包。 tar xzf genlib2-$RELEASE.tar.gz tar xzf dgap-$RELEASE.tar.gz 2、编译genlib2 库: cd genlib2-$RELEASE make 3、重命名这个目录为genlib2 或者创建一个符号链接: cd .. ln -s genlib2-$RELEASE genlib2 4、编译DGAP 程序 cd .. ln -s genlib2-$RELEASE genlib2 为了安装DGAP 只需要从主目录中复制程序dgap 到/usr/bin or /usr/local/bin.。 Sparse Bundle Adhustment 1.5使用指南 Sparse bundle adjustment即稀疏集束调整，现在广泛应用于计算机视觉领域，基本成为最后优化的 标准 excel标准偏差excel标准偏差函数exl标准差函数国标检验抽样标准表免费下载红头文件格式标准下载 算法，就是在已经得到的初始摄像机参数和三维点数据的基础针对投影误差进行优化，得到使得均方投影误差最小意义下的 Motion和Structure。其算法的核心是利用Levenberg-Marquardt算法，由于视觉中问题的特殊性，造成矩阵稀疏，从而针对此特性进行求解。这里介绍的是使用比较多的一个工具函数库，由希腊学者Manolis I.A.Lourakis 和Antonis A. Argyros开发，网址为:。这里我主要翻译下其中接口函数的使用说明: 需要使用的结构体是struct sba_crsm， 定义如下: struct sba_crsm{ int nr,nc; // 稀疏矩阵的行和列 int nnz; //非0元素数量 int *val; //非0元素存储空间，大小为nnz int *colidx //非零元素的列下标，大小为nnz int *rowptr; //val中开始一行的位置，大小为:nr+1, rowptr[nr] = nnz } 本结构体用于保存稀疏矩阵的元素 (根据其定义可以看出，其存储过程其实对于稀疏矩阵按行扫描，得 到非零元素出现的行和列，然后用colidx存储列下标，rowptr存储数据中 新的一行出现的位置，非零结果存在val中。) Sparse Bundle Adjustment通过函数sba_motstr_levmar_x()执行，函数 原型如下: int sba_motstr_levmar_x(const int n, const int m, const int mcon, char * vmask, double *p, const int cnp, const int pnp, double *x, const int mnp, void(*func)(double *p, struct sba_crsm * idxij, int *wk1, int *wk2, double *hx, void *adata), void(*fjac)(double *p, struct sba_crsm * idxij, int *wk1, int *wk2, double *hx, void *adata), void *adata, int itmax, int verbose, double opts[3], double info[10]; 函 数运行成功终止则返回达到最小需要的迭代次数，否则输入-1。当前版本 中假设观测向量的协方差矩阵为单位矩阵。函数的核心是通过LU分解求 解线性系统，使用LAPACK实现。下面是用I,O分别表示输入和输出: n: 3D点的数量。(I) m:摄像机的数量(图像)(I) mcon:从第一个摄像机开始不需要修正参数的摄像机的数量。当我们 将世界坐标系与第一个摄像机的坐标对齐的时候，第一个摄像机矩阵保持 为[I 0] (I) vmask:点的可见性 标志 禁止坐卧标志下载饮用水保护区标志下载桥隧标志图下载上坡路安全标志下载地理标志专用标志下载 :mask[(i-1)*m + j -1] = 1 如果点i在图像j中可见，否则为0.需要注意的是点和摄像机的表示从1，2开始而矩阵下标则是从0开始。Vamsk的大小为n*m(I); p:作为输入，初始参数向量P0 = (a1,a2,?am, b1,b2,?bn)。作为输出，估计的最小值，其大小为m*cnp+n*pnp (I/O); cnp: 定义一个摄像机的参数数量。例如，欧氏摄像机参数化使用6个参数(3个旋转参数+3个平移参数)。如果使用四元数表述旋转则参数的数量增加到7(4+3)。完整的透视摄像机使用11个参数表示，如果包括全局的尺度因子则参数为12个 (I) pnp:定义一个3D 点需要的参数数量，对于欧氏几何为3，对于投影空间为4 (I) x:观测向量X包括所有图像投影，顺序为 (x11,x12,?,x1m,?.xnm).对于不可见的点，对应的xij消失。最大为n*m*mnp.(I) mnp:定义一个图像点参数的数量(一般为2) (I) func:计算估计的参数向量的函数。 (I) fjac:在jac中评估<!--[if !vml]--><!--[endif]-->处的稀疏雅可比行列式<!--[if !vml]--><!--[endif]-->(I) adata:指向可能的额外数据的指针，传递给func，fjac。主要是为了避免对于全局变量的直接应用(I) itmax:Levenberg-Marquardt迭代的最大次数。 (I) verbose:冗长程度。0表示冷静的操作，大的值对应着增加冗长级别。(I) opts:Levenberg-Marquardt算法中最小化参数选项， <!--[if !vml]--> <!--[endif]--> ，分别对应着初始衰减项的尺度因子和结束的终止门限(I); info:关于最小化输出的信息，如果不需要可以设置为NULL。(O) info[0]:<!--[if !vml]--><!--[endif]-->初始参数估计的误差。主要到info[0]除以所有的图像点观测的数量对应着初始均方投影误差; info[1-4]:(<!--[if !vml]--><!--[endif]-->，<!--[if !vml]--> <!--[endif]-->，<!--[endif]-->，<!--[if !vml]--> <!--[if !vml]--> <!--[if !vml]--><!--[endif]-->)，全部是在<!--[endif]-->计算得到的。类似于info[0]，info[1]除以图像观测点的数量得到最终的均方投影误差; info[5]:总的迭代次数; info[6]:结束的原因: 1.过小的<!--[if !vml]--><!--[endif]--> 2 过小的<!--[if !vml]--> 3 迭代次数达到itmax <!--[endif]--> 4 增强法方程矩阵奇异，最小化过程应该从当前解重新开始采用增加的衰减项; Info[7]:function评价的总的次数; Info[8]:fjac评价的总的次数; Info[9]:增强法方程求解总的次数。这个参数一般比总的迭代次数要大，由于在一次LM迭代中，可能需要多个参数进行尝试，每一个都需要对应的增强法方程的求解。 另外Sba还提供了两个函数sba_mot_levmar_x()和sba_str_levmat_x()，分别只对相机参数和结构参数优化投影误差。 个别公式无法显示。 关于sba(sparse bundle adjustment)的30个常见问题(1-14) sba FAQ Q1 -- 什么是sba, sba是一个C/C++软件包对广义稀疏光束平差，在GNU公共许可证下分发。sba是通用的，提供关于定义涉及光束法平差的图像投影的参数选择和函数关系增强的灵活性。 Q2 -- 什么是光束法平差, 假设给定一系列图像中观测到的一组对应点集相应的三维坐标的初始估计，以及关于每张图像的viewing参数的初始估计。光束法平差(BA)是一个大的最优化的问题，包括同时精化三维结构和viewing参数(即相机姿态和可能的本征校准和径向畸变)，为了获得一个在特定的假设下最优化的重建，考虑与观测的图像特征有关的噪声:如果图像误差满足均值为零的正态分布，那么BA是最大似然法估计。它的名字“bundles”(光束)源于每个三维特征聚焦于每个相机的光学中心，这些光学中心相对于结构和viewing参数进行最优化的调整。sba使用Levenberg-Marquardt非线性最小二乘算法的常规实现来解决与BA相联系的稀疏的大规模的优化问题。 Q3 -- “稀疏”是什么意思, 由于不同的三维点和相机之间没有相互影响，BA过程中必须解决的线性系统(即法方程)包含许多零并且形成了一个稀疏的块结构特征。这种 结构可以通过避免存储和处理零元素来利用，从而获得可观的计算效益。 Q4 -- 为什么通用的优化代码不能用于实现BA, BA涉及大规模最小化问题的解决，一般涉及上千个变量。优化算法通过迭代的进行函数的线性化，在当前估计的周围进行最小化并且获得线性系统(它的解确定了当前估计的一个增量)。大多数优化代码(如，minpack, nl2sol/n2g 等)假设这些线性系统是稠密的，即包含很多非零元素。已知的几个这种系统的计算复杂度是O(n**3)。调用BA涉及到大规模最小化问题的解(一般涉及上千个变量)，所以很明显大多数通用的优化代码在应用BA时会引起运行效率降低和存储容量增大。 Q5 -- 在哪里可以找到更多的关于光束法平差的信息, 关于BA的在基于视觉的重建中的应用的一个优秀的综述:Triggs 等人的Bundle Adjustment: A Modern Synthesis。更简短的说明，可参考 bundle adjustment article in Wikipedia。 Q6 -- 在哪里可以找到更多关于Levenberg-Marquardt算法的信息, 参考讲稿Methods for Non-Linear Least Squares Problems, by K. Madsen, H.B. Nielsen and O. Tingleff, Technical University of Denmark, 2004. Q7 -- sba 支持哪些类型的光束法平差, Sba 给它的用户对于描述相机和三维结构的参数定义的完全的控制。因此，它可以支持不同的多视重建问题的实例，如投影重建，几何重建，本征相机参数精化，等等。Sba也提供了程序处理前方交会和后方交会问题，其中相机姿态和场景结构分别保持不变。 Q8 -- 谁在使用sba, Sba对于计算机视觉、机器人、基于图像的图形学、摄影测量、测量、制图等领域的研究者和从业者是无价的。它已经在全球许多实验室中使用，并且也是世界范围内以源代码方式提供的GNU GPL许可证而且授权商用的软件。如果你想要知道更多关于sba的用途方面的信息，以下是我们列出的使用了sba的论文列表here。 Q9 -- 使用sba时我需要什么, 为了能够使用所有的sba函数，必须安装LAPACK或者等价的库，检查 的f2c'ed的免费版本。据 报告 软件系统测试报告下载sgs报告如何下载关于路面塌陷情况报告535n,sgs报告怎么下载竣工报告下载 上述站点上的预编译的MSWin库损坏了，需要使用工程文件重新构建。你也可以尝试这些预编译的MSWin LAPACK/BLAS 库。 Q10 -- 我怎么样编译sba, 首先，你必须保证在本地安装了LAPACK。如果你必须安装它，请遵循以下安装说明包括LAPACK的分发。第二步是编译sba本身。tar压缩包包含Unix/Linux下使用gcc的makefile文件以及MSWin下使用Visual Studio 的Makefile.vc。请阅读这些文件中的注释以获得更多信息。基于提供的makefile文件，它可以使用任何符合ANSI的编译器来直接编译。 Q11 -- 为什么我得到一个未确定的外部符号dgesdd的链接错误, 可能是因为你的LAPACK是旧的版本。Dgesdd在LAPACK3.0中引入。如果你不想升级，你可以选择使用稍微慢一点的dgesvd。Sba_lapack.c对于注释中的dgesvd有类似的调用。 Q12 -- 不同的事物参数在哪里有详细的解释, 参考 ICS/FORTH TR-340: The Design and Implementation of a Generic Sparse Bundle Adjustment Software Package Based on the Levenberg-Marquardt Algorithm, by M.I.A. Lourakis and A.A. Argyros, 2004. 源代码中也包含了每个函数的每个参数的注释。 Q13 -- 怎样把sba改编为我自己的BA变体, 如上所述，sba在选择描述相机、三维结构和图像投影的函数关系和 参数时非常灵活。因此，它可以支持许多种类的BA，包括任意投影或者仿 射相机，部分地活着完全的本征标定相机，外方位元素(即姿态)估计， 从特定的三维点、本征标定图像的三维重建，本征标定参数的精化，等等。 为了把sba改编为一个特定的问题，用户必须对相机和三维结构选择一个 合适的参数化方法，然后提供相应的投影函数以及它的函数行列式的实现 代码。包括sba的演示程序更详细的阐述了处理几何BA问题的细节。 Q14 -- 如何计算sba重建的初始点, 由于BA涉及一个迭代最小化问题的解，不同参数的初始估计应该提 供给sba。这种估计定义了一个初始的三维重建并且可以使用任何结构或 者运动估计视觉算法来计算，如 Hartley and Zisserman's 的书里面描 述的方法。 本文全部译自~lourakis/sba/faq.html Q15 -- sba对outlying数据稳健吗, 短的回答:No。 长的回答:假设图像点特征的局部误差时零平均的正态分布，BA是最 大似然法估计。正是由于这个性质，使得BA如此强大，同时精化三维结 构和viewing 参数。然而，注意，前述的假设不包括有gross outliers(即 不匹配的特征)的情况。这些outlying特征应该在BA之前使用几何方法检测和消除。这种方法的更多的细节，在这里。 Q16 -- 怎样避免指定事务来计算投影函数行列式, 通过传入NULL作为函数行列式的对应参数，函数行列式在前向有限差分方法的辅助下进行计算。然而，注意，我们不推荐这种选择:为了获得最大的效率，我们建议提供一个函数解析的评估函数行列式。 Q17 -- 我们如何通过投影函数和函数行列式来传入和使用我们自己的数据, 设想你想传入两个数组，一个是double型的，一个是integer型的。一个快速并且有效的方法是使用全局变量。避免使用全局变量的最简单的方法是声明一个结构体如下 struct mydata{ double dar[XXX]; int iar[YYY]; }; 其中XXX和YYY表示合适的数组大小。然后，定义一个结构体变量: struct mydata data; 然后赋值: data.dar[0]=7.0; data.iar[0]=-17; // etc 然后，调用合适sba程序，传入数据的地址作为参数，如: ret=sba_motstr_levmar(..., proj, projac, (void *)&data, itmax, verbose, opts, info); 你的proj 和 projac 程序需要使用类型转换来解析提供的数据: struct mydata *dptr; dptr=(struct mydata *)adata; // adata is passed as void * // supplied data can now be accessed as dptr->dar[0], etc Q18 -- 我如何检验用户提供的jacobi行列式, 文件sba_chkjac.c包括检查传入BA函数的函数行列式的正确性(即用户提供的投影函数的一致性)函数。这些函数可以直接调用，通过将所有的 sba_XXX_levmar() 和 sba_XXX_levmar_x() 中的itmax的值指定为0。这种情况下，这些函数立即返回0。检查函数输出可能的梯度(即函数行列式的行)到stderr。然而，我们应该注意，这些函数并不是100%安全可靠的，因为他们依赖点的评估，他们可能报告函数行列式是正确的。解决这种问题的方法是试着检查其它点Point(s)处的函数行列式的正确性。参考sba_chkjac.c查看更多细节。同时，记住这些函数是为了检查解析函数行列式，而不是有限差分的数值近似。 Q19 -- 为什么我获得一个LAPACK错误:函数sba_Axb_Chol()中错误的因式分解, 完整的消息是这样的: LAPACK error: the leading minor of order XX is not positive definite,the factorization could not be completed for dpotf2/dpotrf insba_Axb_Chol()。这意味着增广的法方程的基于Cholesky的算子在某次迭代时出错。理论上，待解决的系统矩阵应该是正定的，这样Cholesky是有效 的。然而由于数值误差，并不完全是这样的，并且这时Cholesky也无法因式分解。尽管如此，不会有什么损害，因为当一次迭代出错时，法方程的阻尼增加了，所以在接下去的迭代中它们的矩阵式正定的。换言之，这些错误不会影响收敛。为了阻止sba报错，注释掉sba_levmar.c中的对sba_Axb_Chol()的调用，而去掉sba_Axb_LU()的注释。这样我们就能够使用LU分解来取代Choleskyl来解法方程，并且可以避免警告。然而，注意，LU分解需要更长时间计算，所以以上更改会使sba运行的稍慢一些。 Q20 -- 为什么sba返回错误代码7(用户错误), 这个错误在sba函数探测到用户生成的预测投影的一个或多个无效的(即，NaN 或者Inf)值返回。如果你使用自己的投影函数，那么就要保证他们被正确的编码了。对于投影函数，包括sba的演示程序，一般是由于对特定点的无效的初始估计;然而，它也可能是因为错误的初始相机参数。为了将引起问题的投影打印在屏幕上，重新运行eucsbademo，使verbosity level为2或者更高。这可以通过设置 eucsbademo.c中的verbose变量为2，然后重新编译来达到。 Q21 -- 演示程序解决哪些类型的光束法平差问题, 如Q7中所述，sba在光束法平差的参数化选择方面非常灵活。为了向用户展示几何光束法平差，demo目录包含了educsbademo程序。升级到1.3版本之前，eucsbademo只支持所有图像的相机本征参数相同的并且在光束法平差过程中保持不变的情况。从版本1.3开始，eucsbademo也支持图像的相机本征参数不同的并且在光束法平差过程中变化的情况。 Q22 -- 怎样更改本征参数的个数并在演示程序过程中保持不变, eucsbademo程序支持BA的多种的相机本征参数，每个相机都不同(见Q21)。Eucsbademo支持不同数量的本征参数，并且在BA过程中保持不变。如:歪曲率保持不变，纵横比和歪曲率保持不变，纵横比、歪曲率和主点保持不变，更多细节可以在sba_driver()中的第一行找到。 Q23 -- 演示程序的输出在哪里, 一般情况下演示程序不产生输出;只有关于最小化的一些统计资料会被打印出来。为了在BA之后打印运动和/或结构参数到stdout，请编辑eucsbademo.c中的函数sba_driver()，并且去掉对prnt变量赋值的注释，将prnt赋值为BA_MOTSTRUCT, BA_MOT 或者 BA_STRUCT。 Q24 -- 演示程序中的协定坐标系统是什么, 演示程序假设一个左手相机坐标系统，Z轴与相机光轴方向一致，并且指向场景。另一个左手系是世界坐标系统。请参考这里和这里观看图示。 Q25 -- 在哪里可以找到更多关于四元数的信息, Eucsbademo演示程序包括sba使用四元数来表示旋转。更多关于四元数和旋转的细节可以在Berthold Horn的笔记中找到。 Q26 -- 如何计算图像投影的协方差, 一个图像点的协方差只能够近似的计算;最常见的做法是使用图像深度的空间派生物。更多细节在这里given by Brooks et al. in What value covariance information in estimating vision parameters?, ICCV01, vol. I, pp. 302-308, Vancouver, IEEE Press, 2001. Sba 迭代停止的原因 1 - stopped by small gradient J^T e 2 - stopped by small dp 3 - stopped by itmax 4 - stopped by small relative reduction in ||e||_2 5 - stopped by small ||e||_2 6 - too many attempts to increase damping. 增大mu，重启sba。 7 - stopped by invalid (i.e. NaN or Inf) "func" values. 投影函数或数据错误。 Sba 中的投影函数 Sba 的投影方程为 其中A即为sba中需要的相机标定矩阵(3′3)，在传入的相机标定文件 (calib.txt)中给出;K 为sba中的旋转矩阵(3′3)，由相机文件(cams.txt)给出其对应的四元数;t为sba中的平移量(3′1)，在相机文件(cams.txt)中直接给出。也可以修改sba 的投影函数calcImgProj，它文件imgproj.c 中。 2009-08-19 19:48 sba中非命令行参数的设置 在sba 工程中配置好传入的参数命令行参数为cams.txt pts.txt calib.txt，函数 sba_driver 中的几个变量的设置如下: nconstframes=0;//nconstframe 表示前几幅图像的参数在平差过程中保持不变。这里我们 使所有图像都参与平差，因此设nconstframe 为0。 howto= BA_MOTSTRUCT;//howto 表示平差的方式，只对相机参数(motion)进行平差,只对 激光点坐标(structure)进行平差,还是一起平差,这里我们选择一起平差。 expert=1;//sba 过程提供两个接口，一个是专家版(运算过程复杂，时间长，结果更精 确)，以及普通版。这里我们选择专家版。 analyticjac=1;//使用解析方法计算投影方程的雅可比矩阵，还是近似计算雅可比矩阵, 这里我们选择解析方法，因为更准确。 prnt = BA_MOTSTRUCT;//打印相机参数的估计值，还是激光点参数的估计值,这里我 们选择同时打印。 计算过程如下:将sba中两个重要的宏SBA_INIT_MU 和SBA_STOP_THRESH置为默认值，或输出的图像点误差仍然不够小，或显示的迭代停止的原因代号为4，则减小SBA_INIT_MU，并将迭代的结果作为初值代入重新平差，直到输出的图像点误差足够小并且迭代停止的原理代号为1或2。以下是一组计算实例: 首先，将宏SBA_INIT_MU的值置为1E-03，SBA_STOP_THRESH的值置为1E-18。运行sba，得到如下结果(将平差值写入文件中): Starting BA with fixed intrinsic parameters SBA using 19943 3D pts, 5 frames and 41563 image projections, 59859 variables Method BA_MOTSTRUCT, expert driver, analytic Jacobian, fixed intrinsics, without covariances SBA returned 6 in 6 iter, reason 2, error 827.631 [initial 2940.98], 6/6 func/fjac evals, 6 lin. systems Elapsed time: 1.06 seconds, 1063.00 msecs 图像点误差为827.631，比较大。将宏SBA_INIT_MU的值置为1E-06， 将运算结果作为初值，重新运行sba，得到如下结果: Starting BA with fixed intrinsic parameters LAPACK error: the leading minor of order 27 is not positive definite, the factorization could not be completed for dpotf2/dpotrf in sba_Axb_Chol() LAPACK error: the leading minor of order 27 is not positive definite, the factorization could not be completed for dpotf2/dpotrf in sba_Axb_Chol() LAPACK error: the leading minor of order 27 is not positive definite, the factorization could not be completed for dpotf2/dpotrf in sba_Axb_Chol() LAPACK error: the leading minor of order 27 is not positive definite, the factorization could not be completed for dpotf2/dpotrf in sba_Axb_Chol() LAPACK error: the leading minor of order 27 is not positive definite, the factorization could not be completed for dpotf2/dpotrf in sba_Axb_Chol() LAPACK error: the leading minor of order 27 is not positive definite, the factorization could not be completed for dpotf2/dpotrf in sba_Axb_Chol() LAPACK error: the leading minor of order 27 is not positive definite, the factorization could not be completed for dpotf2/dpotrf in sba_Axb_Chol() LAPACK error: the leading minor of order 27 is not positive definite, the factorization could not be completed for dpotf2/dpotrf in sba_Axb_Chol() SBA using 19943 3D pts, 5 frames and 41563 image projections, 59859 variables Method BA_MOTSTRUCT, expert driver, analytic Jacobian, fixed intrinsics, without covariances SBA returned 1369 in 1369 iter, reason 4, error 0.300721 [initial 827.631], 1381 /1369 func/fjac evals, 1388 lin. systems Elapsed time: 246.28 seconds, 246282.00 msecs 图像点误差为0.300721，已经足够小，但迭代停止的原因代号为4， 所以将平差结果重新传入，重新运行sba，得到如下结果: Starting BA with fixed intrinsic parameters SBA using 19943 3D pts, 5 frames and 41563 image projections, 59859 variables Method BA_MOTSTRUCT, expert driver, analytic Jacobian, fixed intrinsics, without covariances SBA returned 1 in 1 iter, reason 2, error 0.300721 [initial 0.300721], 1/1 func/ fjac evals, 1 lin. systems Elapsed time: 0.19 seconds, 188.00 msecs 此时，图像点误差为0.300721，已经足够小，并且迭代停止的原因代号为2，所以平差值可以作为最后的结果，平差结束。得到三个文件“pts平差值.txt”、“cams平差值.txt”，“3DPoints平差值.txt”。 使用SBA(Sparse Bundle Adjustment)时应传入的的命令行参数 分两种情况: 1、当所有相机的内参数相同时(使用型号相同的相机或者所有相片使用同一相机拍摄)需传入三个命令行参数: 按顺序以空格或tab符间隔 cams.txt pts.txt calib.txt。 cams.txt为相机外参数文件，描述了每张照片的外参数。格式如下: r a b c Xt Yt Zt 其中r a b c 为投影函数旋转矩阵对应的四元数(r为实部，a b c为虚部)。Xt Yt Zt为平移量。 pts.txt为激光点与图像点的匹配文件，格式如下: X Y Z nframes frame0 x0 y0 frame1 x1 y1 ... 文件的每一行如上所示。X Y Z为激光点坐标。nframes为该激光点对应的图像 点的个数。后面是对nframes个图像点的描述。frame0 x0 y0 表示第frame0幅照片上图像坐标为x0 y0的像点，以此类推。 calib.txt由相机规格和相片大小决定的相机标定文件。跟计算机视觉中的相机标定矩阵一致。可以根据相机焦距、像元大小以及相片的宽度和高度写出相机标定文件calib.txt。比如，我们获得的焦距f = 28.1359mm， 像元大小为8u，相片宽度为3000像素，高度为4500像素。那么我们得到的相机标定文件为:(通常不考虑相机坐标轴不平行的情况) 3516.9875 0.0 1500 0.0 3516.9875 2250 0.0 0.0 1.0 2、当相机的内参数不同时，需传入两个命令行参数: camsvarK.txt pts.txt 因为相机的内参数不同，所以没有统一的相机标定文件。 这里，参数pts.txt不变。camsvarK.txt表示了每张照片对应的相机的内外参数。格式如下: fu, u0, v0, ar, s quaternion translation fu表示以像素为单位的焦距大小，u0,v0 表示像主点坐标。 比如，若焦距f = 28.1359mm，像元大小为8u，相片宽度为3000像素，高度为4500像素，则fu = 3516.9875,u0 = 1500,v0 = 2250。 通常设ar=1,s=0。 ************************************************************** SBA version 1.6 By Manolis Lourakis Institute of Computer Science Foundation for Research and Technology - Hellas Heraklion, Crete, Greece ************************************************************** ==================== GENERAL ==================== This is sba, a copylefted C/C++ implementation of generic bundle adjustment based on the sparse Levenberg-Marquardt algorithm. sba can support a wide range of manifestations/parameterizations of the multiple view reconstruction problem such as arbitrary projective cameras, partially or fully intrinsically calibrated cameras, exterior orientation (i.e. pose) estimation from fixed 3D points, etc. sba can be downloaded from ~lourakis/sba sba relies on lapack for solving the augmented normal equations arising in the course of the Levenberg-Marquardt algorithm. if you don't already have lapack, I suggest getting clapack from Directory demo contains eucsbademo, a working example of using sba for Euclidean bundle adjustment. More details regarding sba can be found in ICS/FORTH Technical Report No. 340 entitled "The Design and Implementation of a Generic Sparse Bundle Adjustment Software Package Based on the Levenberg-Marquardt Algorithm", by M.I.A. Lourakis and A.A. Argyros (available from ~lourakis/sba) In case that you use sba in your published work, please include a reference to the above TR: @techreport{lourakis04, author={M.I.A. Lourakis and A.A. Argyros}, title={The Design and Implementation of a Generic Sparse Bundle Adjustment Software Package Based on the Levenberg-Marquardt Algorithm} institution={Institute of Computer Science - FORTH}, address={Heraklion, Crete, Greece}, number={340}, year={2004}, month={Aug.}, note={Available from \verb+~lourakis/sba+} } ==================== FILES ==================== sba_levmar.c: SBA expert driver routines sba_levmar_wrap.c: simple wrappers around the routines in sba_levmar.c sba_lapack.c: LAPACK-based linear system solvers (LU, QR, SVD, Cholesky, Bunch-Kaufman) sba_crsm.c: CRS sparse matrix manipulation routines sba_chkjac.c: routines for verifying the correctness of user-supplied jacobians sba.h: Function prototypes & related data structures demo/*: Euclidean BA demo; see demo/README.txt for more details matlab/*: sba MEX-file interface; see matlab/README.txt for more details utils/*: Various utilities; see utils/README.txt for more details ==================== COMPILING ==================== - On a Linux/Unix system, typing "make" will build both sba and the demo program. - Under Windows and if Visual C is installed & configured for command line use, type "nmake /f Makefile.vc" in a cmd window to build sba and the demo program. In case of trouble, read the comments on top of Makefile.vc Send your comments/bug reports to lourakis@ics.forth.gr