matlab 调用 opengl 画骨架

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


博主写程序一般用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


评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值