博主写程序一般用matlab, 因为最方便, 基本不用写基本函数, 只需要调用写好的库就行了。 然而, matlab并不是万能的, 有些事情matlab就不擅长, 比如画三维图形。
经研究, 可以采用折中的解决方案, 即: 使用opengl画三维图形, 同时离线渲染, 保存到内存中, 然后用matlab读取内存, 得到结果。
这里, 离线渲染是指不显示窗口, 直接渲染到数组里。 所以glut是不能用了, 我采用的是osmesa。 这是一个开源的linux opengl库, 不过也可以在windows下编译。
编译osmesa见http://stackoverflow.com/questions/17871781/building-mesa-for-windows-7-mesa-9-1, 照做基本没问题, 我使用mesa 10.4, vs2010成功编译出了32位和64位。
得到osmesa之后, 可以参考mesa自带的demo, 然后结合matlab的mex函数, 即可在matlab调用opengl。
如示例, 这里我写了一个draw skeleton的程序。
parse。h
#ifndef REPARSE_H_
#define REPARSE_H_
#include <iostream>
using namespace std;
#include <vector>
typedef vector<float> mypoint;// x, y, z
typedef vector<mypoint> myskeleton;
myskeleton readsingle(double *da)
{
int jo = 15;
myskeleton skel;
skel.resize(jo);
for (int j = 0; j < 15; j++)
{
skel[j].resize(3);
skel[j][0] = da[j];
skel[j][1] = da[j + 15];
skel[j][2] = da[j + 30];
}
return skel;
}
#endif
show.cpp
#include <GL/osmesa.h>
#include <gl/GLU.h>
#include <gl/glut.h>
#include <mex.h>
#include <iostream>
using namespace std;
#include "reparse.h"
static myskeleton mk;
static GLfloat width = 640;
static GLfloat height = 480;
// 已知三点坐标,求法向量
float *CalNormal(float *p1, float * p2, float * p3)
{
float a = ( (p2[1]-p1[1])*(p3[2]-p1[2])-(p2[2]-p1[2])*(p3[1]-p1[1]) );
float b = ( (p2[2]-p1[2])*(p3[0]-p1[0])-(p2[0]-p1[0])*(p3[2]-p1[2]) );
float c = ( (p2[0]-p1[0])*(p3[1]-p1[1])-(p2[1]-p1[1])*(p3[0]-p1[0]) );
float vnorm[3] = {a, b, c};
return vnorm;
}
// draw sphere
void drawSphere(GLfloat xx, GLfloat yy, GLfloat zz, GLfloat radius, GLfloat M, GLfloat N)
{
GLfloat PI = 3.1416;
float step_z = PI / M;
float step_xy = 2 * PI / N;
float x[4], y[4], z[4];
float angle_z = 0.0;
float angle_xy = 0.0;
int i = 0, j = 0;
glBegin(GL_QUADS);
for (i = 0; i < M; i++)
{
angle_z = i * step_z;
for (j = 0; j < N; j++)
{
angle_xy = j * step_xy;
x[0] = radius * sin(angle_z) * cos(angle_xy);
y[0] = radius * sin(angle_z) * sin(angle_xy);
z[0] = radius * cos(angle_z);
x[1] = radius * sin(angle_z + step_z) * cos(angle_xy);
y[1] = radius * sin(angle_z + step_z) * sin(angle_xy);
z[1] = radius * cos(angle_z + step_z);
x[2] = radius*sin(angle_z + step_z)*cos(angle_xy + step_xy);
y[2] = radius*sin(angle_z + step_z)*sin(angle_xy + step_xy);
z[2] = radius*cos(angle_z + step_z);
x[3] = radius * sin(angle_z) * cos(angle_xy + step_xy);
y[3] = radius * sin(angle_z) * sin(angle_xy + step_xy);
z[3] = radius * cos(angle_z);
float p0[3] = {x[0], y[0], z[0]};
float p1[3] = {x[1], y[1], z[1]};
float p2[3] = {x[2], y[2], z[2]};
glNormal3fv(CalNormal(p0, p1, p2));
for (int k = 0; k < 4; k++)
{
glVertex3f(xx + x[k], yy + y[k], zz + z[k]);
}
}
}
glEnd();
}
void RenderBone(float x0, float y0, float z0, float x1, float y1, float z1)
{
GLdouble dir_x = x1 - x0;
GLdouble dir_y = y1 - y0;
GLdouble dir_z = z1 - z0;
GLdouble bone_length = sqrt( dir_x*dir_x + dir_y*dir_y + dir_z*dir_z );
// glPushMatrix();
glTranslated( x0, y0, z0 );
double length;
length = sqrt( dir_x*dir_x + dir_y*dir_y + dir_z*dir_z );
if ( length < 0.0001 ) {
dir_x = 0.0; dir_y = 0.0; dir_z = 1.0; length = 1.0;
}
dir_x /= length; dir_y /= length; dir_z /= length;
GLdouble up_x, up_y, up_z;
up_x = 0.0;
up_y = 1.0;
up_z = 0.0;
double side_x, side_y, side_z;
side_x = up_y * dir_z - up_z * dir_y;
side_y = up_z * dir_x - up_x * dir_z;
side_z = up_x * dir_y - up_y * dir_x;
length = sqrt( side_x*side_x + side_y*side_y + side_z*side_z );
if ( length < 0.0001 ) {
side_x = 1.0; side_y = 0.0; side_z = 0.0; length = 1.0;
}
side_x /= length; side_y /= length; side_z /= length;
up_x = dir_y * side_z - dir_z * side_y;
up_y = dir_z * side_x - dir_x * side_z;
up_z = dir_x * side_y - dir_y * side_x;
GLdouble m[16] = { side_x, side_y, side_z, 0.0,
up_x, up_y, up_z, 0.0,
dir_x, dir_y, dir_z, 0.0,
0.0, 0.0, 0.0, 1.0 };
glMultMatrixd( m );
GLdouble radius= 0.01;
glBegin(GL_TRIANGLES);
float d1[3] = {-1*radius, -1*radius, 0*bone_length};
float d2[3] = {1*radius, -1*radius, 0*bone_length};
float d3[3] = {1*radius, 1*radius, 0*bone_length};
float d4[3] = {-1*radius, 1*radius, 0*bone_length};
float u1[3] = {-1*radius, -1*radius, 1*bone_length};
float u2[3] = {1*radius, -1*radius, 1*bone_length};
float u3[3] = {1*radius, 1*radius, 1*bone_length};
float u4[3] = {-1*radius, 1*radius, 1*bone_length};
glNormal3f(0, 0, 1);
glVertex3fv(u1);
glVertex3fv(u2);
glVertex3fv(u3);
glVertex3fv(u1);
glVertex3fv(u3);
glVertex3fv(u4);
glNormal3f(1, 0, 0);
glVertex3fv(u2);
glVertex3fv(d2);
glVertex3fv(d3);
glVertex3fv(u2);
glVertex3fv(d3);
glVertex3fv(u3);
glNormal3f(0, 0, -1);
glVertex3fv(d1);
glVertex3fv(d4);
glVertex3fv(d3);
glVertex3fv(d1);
glVertex3fv(d3);
glVertex3fv(d4);
glNormal3f(-1, 0, 0);
glVertex3fv(u1);
glVertex3fv(u4);
glVertex3fv(d4);
glVertex3fv(u1);
glVertex3fv(d4);
glVertex3fv(d1);
glNormal3f(0, 1, 0);
glVertex3fv(u4);
glVertex3fv(u3);
glVertex3fv(d3);
glVertex3fv(u4);
glVertex3fv(d3);
glVertex3fv(d4);
glNormal3f(0, -1, 0);
glVertex3fv(u1);
glVertex3fv(d1);
glVertex3fv(d2);
glVertex3fv(u1);
glVertex3fv(d2);
glVertex3fv(u2);
glEnd();
GLdouble m2[16] = { side_x, up_x, dir_x, 0.0,
side_y, up_y, dir_y, 0.0,
side_z, up_z, dir_z, 0.0,
0.0, 0.0, 0.0, 1.0 };
glMultMatrixd( m2 );
glTranslated( -x0, -y0, -z0 );
// glPopMatrix();
}
void Display(void)
{
glClear(GL_COLOR_BUFFER_BIT|GL_DEPTH_BUFFER_BIT);
glPushMatrix();
glColor3f(1.0, 0.0, 0.0);
for (int i = 0; i< 15; i++)
{
mypoint mp = mk[i];
glTranslatef(mp[0], mp[1], mp[2]);
// glutSolidSphere(0.06, 10, 10);
drawSphere(0, 0, 0, 0.03, 10, 10);
glTranslatef(-mp[0], -mp[1], -mp[2]);
}
glColor3f(1, 1, 0);
float ind[] = {1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14};
float ind2[] = {0, 1, 2, 3, 1, 5, 6, 1, 8, 9, 10, 8, 12, 13};
// float ind[] = {1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16};
// float ind2[] = {0, 1, 2, 3, 1, 5, 6, 1, 8, 9, 10, 11, 8, 13, 14, 15};
for (int i = 0; i< 14; i++)
{
int in1 = ind[i];
int in2 = ind2[i];
mypoint mp = mk[in1];
mypoint mp2 = mk[in2];
// glVertex3f(mp[0], mp[1], mp[2]);
// glVertex3f(mp2[0], mp2[1], mp2[2]);
if ( i == 0 || i == 7)
{
glColor3f(1, 1, 1);
}
if ( i == 1 || i == 2 ||
i == 3 || i == 8 ||
i == 9 || i == 10 )
{
glColor3f(0, 0.5, 0.7);
}
if ( i == 4 || i == 5 ||
i == 6 || i == 11 ||
i == 12 || i == 13 )
{
glColor3f(0, 1, 0);
}
RenderBone(mp[0], mp[1], mp[2], mp2[0], mp2[1], mp2[2]);
}
glPopMatrix();
}
void SetupRC(void)
{
glClearColor (1, 1, 1, 1);
glClearDepth(1);
glEnable(GL_DEPTH_TEST);
glEnable(GL_NORMALIZE);
glShadeModel(GL_SMOOTH);
GLfloat _ambient[]={1.0,1.0,1.0,1.0};
GLfloat _diffuse[]={1.0,1.0,1.0,1.0};
GLfloat _specular[]={1.0,1.0,1.0,1.0};
GLfloat _position[]={0,0,200,0};
glLightfv(GL_LIGHT0,GL_AMBIENT,_ambient);
glLightfv(GL_LIGHT0,GL_DIFFUSE,_diffuse);
glLightfv(GL_LIGHT0,GL_SPECULAR,_specular);
glLightfv(GL_LIGHT0,GL_POSITION,_position);
glEnable(GL_LIGHT0);
glEnable(GL_LIGHTING);
glColorMaterial(GL_FRONT, GL_DIFFUSE);
glEnable(GL_COLOR_MATERIAL);//启用颜色追踪
// best quality
glHint(GL_PERSPECTIVE_CORRECTION_HINT, GL_NICEST);
}
void render_image()
{
GLfloat nRange = 1.0f;
// Reset coordinate system
glMatrixMode(GL_PROJECTION);
glLoadIdentity();
// Establish clipping volume (left, right, bottom, top, near, far)
if (width <= height)
glOrtho(-nRange, nRange, -nRange * height / width, nRange * height / width, -nRange*10, nRange*10);
else
glOrtho(-nRange * width / height, nRange * width / height, -nRange, nRange, -nRange*10, nRange*10);
glMatrixMode(GL_MODELVIEW);
glLoadIdentity();
// lights
SetupRC();
// show
Display();
/* This is very important!!!
* Make sure buffered commands are finished!!!
*/
glFinish();
}
void get_pic(double *da, double *test1, double *test2, double *test3)
{
mk = readsingle(da);
for (int i = 0; i< mk.size(); i++)
{
mexPrintf("%f %f %f\n", mk[i][0], mk[i][1], mk[i][2]);
}
OSMesaContext ctx;
GLubyte *buffer;
/* Create an RGBA-mode context */
#if OSMESA_MAJOR_VERSION * 100 + OSMESA_MINOR_VERSION >= 305
/* specify Z, stencil, accum sizes */
ctx = OSMesaCreateContextExt(OSMESA_RGBA, 16, 0, 0, NULL);
#else
ctx = OSMesaCreateContext(OSMESA_RGBA, NULL);
#endif
if (!ctx)
{
printf("OSMesaCreateContext failed!\n");
return;
}
/* Allocate the image buffer */
buffer = (GLubyte*)malloc(width * height * 4 * sizeof(GLubyte));
if (!buffer)
{
printf("Alloc image buffer failed!\n");
return;
}
/* Bind the buffer to the context and make it current */
if (!OSMesaMakeCurrent(ctx, buffer, GL_UNSIGNED_BYTE, width, height))
{
printf("OSMesaMakeCurrent failed!\n");
return;
}
{
int z, s, a;
glGetIntegerv(GL_DEPTH_BITS, &z);
glGetIntegerv(GL_STENCIL_BITS, &s);
glGetIntegerv(GL_ACCUM_RED_BITS, &a);
// printf("Depth=%d Stencil=%d Accum=%d\n", z, s, a);
}
render_image();
int y, x, i;
const GLubyte *ptr = buffer;
for (y = height - 1; y >= 0; y--)
{
for (x = 0; x < width; x++)
{
i = (y * width + x) * 4;
int index = x * height + height - 1 - y;
test1[index] = ptr[i];
test2[index] = ptr[i + 1];
test3[index] = ptr[i + 2];
}
}
mexPrintf("all done\n");
/* free the image buffer */
free(buffer);
/* destroy the context */
OSMesaDestroyContext(ctx);
return;
}
void mexFunction(int nlhs, mxArray *plhs[], int nrhs, const mxArray *prhs[])
{
/* Check for proper number of arguments */
// [r g b] = mesa_skel( skel )
// 15x3 skeletons
if (nrhs != 1)
{
mexErrMsgIdAndTxt("MATLAB:mexcpp:nargin",
"MEXCPP requires 1 input arguments.");
}
// r, g, b
else if (nlhs != 3)
{
mexErrMsgIdAndTxt("MATLAB:mexcpp:nargout",
"MEXCPP requires 3 output argument.");
}
// skel
double *data;
int M, N;
data = mxGetPr(prhs[0]);
M = mxGetM(prhs[0]);
N = mxGetN(prhs[0]);
mexPrintf("rows : %d\tcols : %d\n", M, N);
if (N != 3 || M != 15)
{
mexErrMsgIdAndTxt("MATLAB:mexcpp:nargin",
"MEXCPP requires 15x3 skeletons");
}
// window
width = 640;
height = 480;
// output
plhs[0] = mxCreateDoubleMatrix(height, width, mxREAL);
plhs[1] = mxCreateDoubleMatrix(height, width, mxREAL);
plhs[2] = mxCreateDoubleMatrix(height, width, mxREAL);
double *outData1 = (double *)mxGetPr(plhs[0]);
double *outData2 = (double *)mxGetPr(plhs[1]);
double *outData3 = (double *)mxGetPr(plhs[2]);
// run
get_pic(data, outData1, outData2, outData3);
return;
}
然后再matlab调用mex, 记得链接应有的库:
mex opengl_show.cpp ...
-I"D:\Program Files\OSMesa\include" ...
-L"D:\Program Files\OSMesa\x86" ...
-losmesa
1, 这里骨架我使用的是画一个盒子来模拟。 其实我本想用glu去画cylinder, 但是画出之后在osmesa中显示不出来, 不知道为什么。
2, 如果用了glu, 记得在32位下做, 不然链接会有问题。
over

博主分享了如何在MATLAB中利用OpenGL进行三维图形绘制,特别是如何通过离线渲染将结果保存到内存,然后在MATLAB中读取。由于MATLAB不擅长三维图形,博主选择了osmesa作为OpenGL库,提供了编译osmesa的链接,并说明了在MATLAB中使用mex函数调用OpenGL绘制骨架的过程。示例程序包含parse.h和show.cpp两个文件,但尝试使用glu画圆柱体时遇到了显示问题。

1090

被折叠的 条评论
为什么被折叠?



