使用OpenGL对单幅照片进行三维重建

用过Google Sketchup的大概都知道Sketchup有个功能换做照片匹配，photo match，就是根据单幅照片重建3d模型，如下图所示：

实际上即通过用户在图片上绘制两组平行线来确定两个灭点，再通过灭点来计算相关参数。所以隐含的条件是选用的照片必须是两点透视的。另外还包含了对相机模型的一些条件设定，比如才用针孔模型，相机的投影中心在照片平面上的投影是照片的中心，图像垂直方向和水平方向的比例尺相同，等等。所以你会看到在一些教程里会说需要使用原始照片，从别的照片里抠出来的图片，或者变形的图片不好使。

搜了一下论文，发现国内好几篇论文都引用了“Using Vanishing Point for Camera Calibration and Coarse 3D Reconstruction from a Single Image”这篇文章，于是找来看了看（springerlinker）。大致的思路是先计算灭点坐标，然后求相机焦距，然后求世界坐标系到相机坐标系的旋转矩阵，最后求平移矩阵。推导比较详细，像我这样完全没学过射影几何的也大概能看懂，个别地方看不明白之后默认是使用了定理。

因为只是个简单的Demo，所以下图不是我程序的运行截图，是我绘制出三色轴后手工叠加到照片上对应位置，只是为了验证一下角度是正确的。

程序开始部分的数据处理完全按照前述论文中的过程走，只是对图像坐标做了点处理，把原点从左上角移到了中心位置，这样在做垂直线的时候比较简单。最后得到了一个旋转矩阵Mo-c，按照论文的公式，Mo-c点乘向量U/Ro＝向量U/Rc，即对确定方向的一个向量，其在世界坐标系里的坐标经过Mo-c矩阵变换后就得到了其在相机坐标系里的坐标。比如世界坐标系里的x轴上的单位向量经过Mo-c变换后得到了相机坐标系里的一个向量x'，那么意味着相机坐标系的x轴上的单位向量经过Mo-c变换后也是向量x'，因为不管哪个坐标系的x轴上的单位向量都是(1,0,0)。所以对应的OpenGL变换如下：

gl.glRotatef(180, 1, 0, 0); // 旋转，得到ijk
gl.glTranslatef(s*OA[0], s*OA[1], s*OA[2]); // 平移
gl.glMultMatrixf(moc); // 再旋转，得到uvw

这样得到的坐标轴xyz（分别为红色、绿色、蓝色）跟论文里的ijk和uvw是一一对应的。

 

最后一个步骤，就是利用照片高度和前面求出的相机焦距来求 OpenGL对应的视场角：

 

 

float theta = (float) (Math.atan(height/2/f) * 180 / Math.PI * 2); 

glu.gluPerspective(theta,(float)width/(float)height,0.1f,100.0f);

到此，假如JOGL使用的canvas与照片同尺寸，则绘制出来的三个彩轴会很好地贴合到你指定的位置，就如同前面合成效果图里示意的那样。

另外说一下JOGL，国内网页介绍这个的貌似都比较古老了，给的链接都不能用，最好是从它主页右边进入Tutorial，有基础的3个示例，包括配置环境和一个简单的demo。

 

 

以下是代码

import java.awt.Frame; import java.awt.event.WindowAdapter; import java.awt.event.WindowEvent; import java.nio.FloatBuffer; import javax.media.opengl.*; import javax.media.opengl.awt.GLCanvas; import javax.media.opengl.glu.GLU; import com.jogamp.opengl.util.Animator; import com.jogamp.opengl.util.FPSAnimator; import com.jogamp.opengl.util.gl2.GLUT; public class JOGLTest implements GLEventListener { private static float f; private static FloatBuffer moc; private static float[] OA; public static void main(String[] args) { Point p0 = new Point(310,529); Point p1 = new Point(353, 166); Point p2 = new Point(593,234); Point p3 = new Point(339,526); Point p4 = new Point(536,499); Point p5 = new Point(318,206); Point p6 = new Point(128,256); Point p7 = new Point(311,530); Point p8 = new Point(86,489); // 两对平行线 Line line1 = new Line(p1,p2); Line line2 = new Line(p3,p4); Line line3 = new Line(p5,p6); Line line4 = new Line(p7,p8); // 两个灭点 Point Fu = new Point(); Point Fv = new Point(); LineCrossPoint(line1, line2, Fu); LineCrossPoint(line3, line4, Fv); // 计算焦距 Point P = new Point(450,324); Point Puv = new Point(); Line Horizon = new Line(Fu, Fv); Line PtoPuv = new Line(-1/Horizon.a, 0); LineCrossPoint(Horizon, PtoPuv, Puv); float PPuv = Distance(P, Puv); float FvPuv = Distance(Fv, Puv); float PuvFu = Distance(Puv, Fu); f = (float) Math.sqrt(FvPuv*PuvFu-PPuv*PPuv); // 计算旋转矩阵 float s1 = (float) Math.sqrt(Fu.x*Fu.x+Fu.y*Fu.y+f*f); float s2 = (float) Math.sqrt(Fv.x*Fv.x+Fv.y*Fv.y+f*f); float uu[] = new float[]{Fu.x/s1, Fu.y/s1, f/s1}; float vv[] = new float[]{Fv.x/s2, Fv.y/s2, f/s2}; float ww[] = new float[3]; ChaCheng(uu,vv,ww); float _moc[] = new float[]{ uu[0], uu[1], uu[2], 0, vv[0], vv[1], vv[2], 0, ww[0], ww[1], ww[2], 0, 0, 0, 0, 1 }; moc = FloatBuffer.wrap(_moc); // 计算平移向量 float s3 = (float) Math.sqrt(p0.x*p0.x+p0.y*p0.y+f*f); OA = new float[]{p0.x/s3,p0.y/s3,f/s3}; // 以下是JOGL配置 GLProfile.initSingleton(); GLProfile glp = GLProfile.getDefault(); GLCapabilities caps = new GLCapabilities(glp); GLCanvas canvas = new GLCanvas(caps); Frame frame = new Frame("AWT Window Test"); frame.setSize(300, 300); frame.add(canvas); frame.setVisible(true); // by default, an AWT Frame doesn't do anything when you click // the close button; this bit of code will terminate the program when // the window is asked to close frame.addWindowListener(new WindowAdapter() { public void windowClosing(WindowEvent e) { System.exit(0); } }); canvas.addGLEventListener(new JOGLTest()); Animator animator = new FPSAnimator(canvas, 60); animator.add(canvas); animator.start(); } // 求两直线焦点 private static void LineCrossPoint(Line line1, Line line2, Point c) { c.x = (int) ((line1.b-line2.b)/(line2.a-line1.a)); c.y = (int) ((line2.a*line1.b-line1.a*line2.b)/(line2.a-line1.a)); } // 求两点平面距离 private static float Distance(Point a, Point b) { return (float) Math.sqrt(Math.pow(a.x-b.x, 2) + Math.pow(a.y-b.y, 2)); } // 叉乘 private static void ChaCheng(float[] a, float[] b, float[] result) { result[0] = a[1]*b[2]-a[2]*b[1]; result[1] = a[2]*b[0]-a[0]*b[2]; result[2] = a[0]*b[1]-a[1]*b[0]; } @Override public void display(GLAutoDrawable arg0) { // TODO Auto-generated method stub GL2 gl = arg0.getGL().getGL2(); GLU glu = GLU.createGLU(gl); GLUT glut = new GLUT(); gl.glClear(GL.GL_COLOR_BUFFER_BIT | GL.GL_DEPTH_BUFFER_BIT); gl.glLoadIdentity(); float s = 20; // 比例系数 gl.glRotatef(180, 1, 0, 0); // 旋转，得到ijk gl.glTranslatef(s*OA[0], s*OA[1], s*OA[2]); // 平移 gl.glMultMatrixf(moc); // 再旋转，得到uvw // 画三轴 gl.glBegin(GL.GL_LINES); gl.glColor3f(1.0f, 0, 0); gl.glVertex3f(0, 0, 0); gl.glVertex3f(3.0f, 0, 0); gl.glEnd(); gl.glBegin(GL.GL_LINES); gl.glColor3f(0, 1.0f, 0); gl.glVertex3f(0, 0, 0); gl.glVertex3f(0, 3.0f, 0); gl.glEnd(); gl.glBegin(GL.GL_LINES); gl.glColor3f(0, 0, 1.0f); gl.glVertex3f(0, 0, 0); gl.glVertex3f(0, 0, 3.0f); gl.glEnd(); } @Override public void dispose(GLAutoDrawable arg0) { // TODO Auto-generated method stub } @Override public void init(GLAutoDrawable arg0) { // TODO Auto-generated method stub GL2 gl = arg0.getGL().getGL2(); gl.glShadeModel(GL2ES1.GL_SMOOTH); // Enables Smooth Shading gl.glClearColor(0.0f, 0.0f, 0.0f, 0.0f); // Black Background gl.glClearDepth(1.0f); // Depth Buffer Setup gl.glEnable(GL.GL_DEPTH_TEST); // Enables Depth Testing gl.glDepthFunc(GL.GL_LEQUAL); // The Type Of Depth Test To Do gl.glHint(GL2ES1.GL_PERSPECTIVE_CORRECTION_HINT, GL.GL_NICEST); // Really Nice Perspective Calculations } @Override public void reshape(GLAutoDrawable arg0, int arg1, int arg2, int width, int height) { // TODO Auto-generated method stub GL2 gl = arg0.getGL().getGL2(); GLU glu = GLU.createGLU(gl); if (height==0) // Prevent A Divide By Zero By { height=1; // Making Height Equal One } gl.glViewport(0, 0, width, height); // Reset The Current Viewport gl.glMatrixMode(GL2ES1.GL_PROJECTION); // Select The Projection Matrix gl.glLoadIdentity(); // Reset The Projection Matrix glu.gluPerspective(45.0f,(float)width/(float)height,0.1f,100.0f); // Calculate The Aspect Ratio Of The Window gl.glMatrixMode(GL2ES1.GL_MODELVIEW); // Select The Modelview Matrix gl.glLoadIdentity(); // Reset The Modelview Matrix } } class Point { public int x, y; public Point() { x=y=0; } public Point(int tx, int ty) { x = tx-450; y = ty-324; } } class Line { public float a, b; public Line(float ta, float tb) { a = ta; b = tb; } public Line(Point p1, Point p2) { a = ((float)(p2.y-p1.y))/(p2.x-p1.x); b = ((float)(p1.y*p2.x-p2.y*p1.x))/(p2.x-p1.x); } }