[討論] 關於mex file的問題...
由於我的程式需要跑大量for loop（六層），
所以我用mex改到c下編寫，
但是c我只有一些基礎（我是物化背景），
加上Matlab裡的mex file又有特別的編寫要求，
在請教很多人以及查了很多資料後有了一些進展，
但是目前碰到的bug很奇怪，
我在for loop裡面做矩陣運算，主要是矩陣的旋轉運算，
隨著我for loop給的角度不同，
而產生隨角度改變的旋轉矩陣，
然而我發現我給的input角度一樣，
output出來的旋轉array不一樣，
像是這樣:
for (i1 = 0; i1 < ntheta; i1++)
{
for (i2 = 0; i2 < nphi; i2++)
{
rotation(R7, phi[i2], theta[i1], 0);
printf("R7 = \n %f %f \n %f %f %f \n %f %f %f \n %f %f %f \n \n"...與下面同行
,phi[i2], theta[i1], R7[0],R7[1],R7[2],R7[3],R7[4],R7[5],R7[6],R7[7],R7[8]);
rotation(R8, phi[i2], theta[i1], 0);
printf("R8 = \n %f %f \n %f %f %f \n %f %f %f \n %f %f %f \n \n"...與下面同行
,phi[i2], theta[i1], R8[0],R8[1],R8[2],R8[3],R8[4],R8[5],R8[6],R8[7],R8[8]);
}
}
print出來的結果如下:
R7 =
0.000000 3.141593 <---- 這是給的input角度
-1.000000 0.000000 -0.000000
0.000000 1.000000 0.000000 <------這是旋轉矩陣
0.000000 0.000000 -1.000000
R8 =
0.000000 3.141593 <-----同樣input
-0.540302 0.841471 -0.000000
0.841471 0.540302 0.000000 <-----output旋轉矩陣不同
0.000000 -0.000000 -1.000000
R7 =
15.707963 3.141593
0.857553 -0.514395 0.000000
-0.514395 -0.857553 0.000000
0.000000 -0.000000 -1.000000
R8 =
15.707963 3.141593
0.654290 0.756244 0.000000
0.756244 -0.654290 0.000000
0.000000 0.000000 -1.000000
R7 =
0.000000 0.000000
0.793480 0.608596 0.000000
-0.608596 0.793480 0.000000
0.000000 0.000000 1.000000
R8 =
0.000000 0.000000
0.701369 0.712799 0.000000
-0.712799 0.701369 0.000000
0.000000 0.000000 1.000000
R7 =
15.707963 0.000000
-0.763960 -0.645264 0.000000
0.645264 -0.763960 0.000000
0.000000 0.000000 1.000000
R8 =
15.707963 0.000000
-0.722102 0.691786 0.000000
-0.691786 -0.722102 0.000000
0.000000 0.000000 1.000000
我的旋轉矩陣函式如下:
#include <stdio.h>
#include <stdio.h>
#include "mex.h"
#include <math.h>
void rotation(double *output, double a, double b, double c);
void rotation(double *output, double a, double b, double c)
{
double xx[9], yy[9], zz[9];
int i, j, k;
for (i = 0; i<9; i++)
{
xx[i] = 0;
yy[i] = 0;
zz[i] = 0;
}
xx[0] = cos(a); xx[1] = sin(a); xx[2] = 0; xx[3] = -sin(a); xx[4] = cos(a); xx[5] = 0; xx[6] = 0; xx[7] = 0; xx[8] = 1;
yy[0] = cos(b); yy[1] = 0; yy[2] = -sin(b); yy[3] = 0; yy[4] = 1; yy[5] = 0; yy[6] = sin(b); yy[7] = 0; yy[8] = cos(b);
zz[0] = cos(c); zz[1] = sin(c); zz[2] = 0; zz[3] = -sin(c); zz[4] = cos(c); zz[5] = 0; zz[6] = 0; zz[7] = 0; zz[8] = 1;
double xy[9], xyz[9];
int size = 3;
for (i=0; i<size; i++)
{
for (j = 0; j<size; j++)
{
xy[size*i+j] = 0;
xyz[size*i+j] = 0;
output[size*i+j] = 0;
}
}
for ( i = 0; i<size; i++)
{
for ( j = 0; j<size; j++)
{
for ( k = 0; k<size; k++)
{
xy[size*i+j] = xy[size*i+j] + xx[size*i+k]*yy[size*k+j];
}
}
}
for ( i = 0; i<size; i++)
{
for ( j = 0; j<size; j++)
{
for ( k = 0; k<size; k++)
{
xyz[size*i+j] = xyz[size*i+j] + xy[size*i+k]*zz[size*k+j];
output[size*i+j] = xyz[size*i+j];
}
}
}
}
落落長一大篇，因為這個程式對我很重要，
如有不足，請再推文告知，謝謝
--
※ 發信站: 批踢踢實業坊(ptt.cc)
◆ From: 140.114.44.41
※ anjackie:轉錄至看板 C_and_CPP 12/20 17:51
→
12/21 09:34, , 1F
12/21 09:34, 1F
→
12/21 09:36, , 2F
12/21 09:36, 2F
我原始碼蠻大的，這只是一開始的旋轉，想要看的測試碼是指全部嗎？
呼叫rotation時用&R7是指把
rotation(R7, phi[i2], theta[i1], 0);
改成rotation(&R7, phi[i2], theta[i1], 0);嗎？
※ 編輯: anjackie 來自: 140.114.44.41 (12/21 11:23)
→
12/21 12:40, , 3F
12/21 12:40, 3F
→
12/21 12:41, , 4F
12/21 12:41, 4F
→
12/21 12:41, , 5F
12/21 12:41, 5F
→
12/21 12:42, , 6F
12/21 12:42, 6F