スポンサーリンク

| キーワード:

gluProjectionを自前実装する

行列の積を自前実装できたのでgluProjectionを自前実装する。

数式

公式の式をそのまま使うとなぜか公式の結果と一致しない。

https://www.khronos.org/registry/OpenGL-Refpages/gl2.1/xhtml/gluProject.xml

実装例を見てみるとどうも式が若干違うのでそっちに合わせる。

http://collagefactory.blogspot.com/2010/03/gluproject-source-code.html

コード


  //誤差の定義
  template<typename T>struct real_error;
  template<>struct real_error<double> { static constexpr double value = 1e-15; };
  template<>struct real_error<float> { static constexpr float value = (float)1e-7; };
  //! @brief 誤差の閾値未満かをチェック
  //! @param [in] value 検査する値
  //! @retval true 値の絶対値は閾値より小さい
  //! @retval false 値の絶対値は閾値より大きい
  template<typename T>
  inline bool is_error(const T value) {
    return
      abs(value) < real_error<T>::value;
  }

  //! @brief IxJ行列の i,j の要素にアクセスする
  //! @param [in] i 行番号
  //! @param [in] j 列番号
  //! @param [in] I 行数
  inline int matijI(const int i, const int j, const int I) {
    return i + I * j;
  }

  //! @brief IxJ行列 × JxK行列 の乗算関数
  //https://suzulang.com/cpp-multmatrix44/
  //! @param [out] C Cik 計算結果格納先
  //! @param [in] A Aijの行列
  //! @param [in] B Bjkの行列
  //! @param [in] I Aの行列の行数
  //! @param [in] J Aの行列の列数
  //! @param [in] K Bの行列の列数
  //! @return C
  template<
    typename real_tC,
    typename real_tA,
    typename real_tB>
  real_tC* mult_matrixIJK(real_tC* C, const real_tA* A, const real_tB* B, const int I, const int J, const int K) {

    for (int i = 0; i < I; i++)
      for (int k = 0; k < K; k++) {
        C[matijI(i, k, I)] = 0.0;
        for (int j = 0; j < J; j++) {
          C[matijI(i, k, I)] += A[matijI(i, j, I)] * B[matijI(j, k, J)];
        }
      }

    return C;
  }

  template<
    typename real_tC,
    typename real_tA,
    typename real_tB>
  real_tC* mult_matrix_vec41(real_tC* C, const real_tA* A, const real_tB* B) {
    return mult_matrixIJK(C, A, B, 4, 4, 1);
  }
  // gluProject 
  // https://www.khronos.org/registry/OpenGL-Refpages/gl2.1/xhtml/gluProject.xml
  // http://collagefactory.blogspot.com/2010/03/gluproject-source-code.html

  //! @brief gluProjectの自前実装版
  //! @param [in] wldX 空間座標のX
  //! @param [in] wldY 空間座標のY
  //! @param [in] wldZ 空間座標のZ
  //! @param [in] model モデルビュー行列
  //! @param [in] proj プロジェクション行列
  //! @param [in] ビューポート x,y,width,height
  //! @param [out] pixX 画面上のX座標
  //! @param [out] pixY 画面上のY座標
  //! @param [out] pixZ 画面上のZ座標 ※計算上は存在する
  //! @retval true 計算に成功
  //! @retval false 計算に失敗
  template<typename real_t>
  inline bool myProject(
    real_t wldX,
    real_t wldY,
    real_t wldZ,
    const real_t* model,
    const real_t* proj,
    const int* view,
    real_t* pixX,
    real_t* pixY,
    real_t* pixZ
  ) {

    real_t world[4] = { wldX,wldY,wldZ,1.0 };
    real_t tmp[4];
    real_t vdd[4];

    mult_matrix_vec41(tmp, model, world); // tmp = M×v
    mult_matrix_vec41(
      vdd, 
      proj, tmp
    ); // v'' = P×tmp
  
    if (is_error(vdd[3]) == true)
      return false;

    vdd[0] /= vdd[3];
    vdd[1] /= vdd[3];
    vdd[2] /= vdd[3]; //v''' = v'' / v''[3]

    *pixX = view[0] + view[2] * (vdd[0] + 1) / 2.0;
    *pixY = view[1] + view[3] * (vdd[1] + 1) / 2.0;
    *pixZ = (vdd[2] + 1) / 2.0;

    return true;

  }

動作チェック

gluProjectionと同じパラメータを入力し同じ結果が出るかを確認する。

テストコード

int WIDTH = 300;
int HEIGHT = 400;

void disp(void) {

  glViewport(0, 0, WIDTH, HEIGHT);


  glMatrixMode(GL_PROJECTION);
glPushMatrix();
 glLoadIdentity(); gluPerspective(45, 0.5, 0.1, 2.5); glMatrixMode(GL_MODELVIEW); glPushMatrix();
glLoadIdentity();
glTranslated(-1, 0, -2); glRotated(1.5, 1.0, 1.0, 1.0);
  ////作成した変換行列を取得する
  GLdouble  modeld[16];
  glGetDoublev(GL_MODELVIEW_MATRIX, modeld);

  GLdouble  projd[16];
  glGetDoublev(GL_PROJECTION_MATRIX, projd);

  GLint view[4];
  glGetIntegerv(GL_VIEWPORT, view);
  // 三次元座標(入力)
  double 
    objx = (rand() % 1000) / 1000.0, 
    objy = (rand() % 1000) / 1000.0,
    objz = (rand() % 1000) / 1000.0;

  // 二次元ピクセル座標(出力)
  double winx, winy, winz;

  gluProject(
    objx, objy, objz,
    modeld, projd, view,
    &winx, &winy, &winz);

  printf("glu %lf %lf %lf\n", winx, winy, winz);
  myProject(objx, objy, objz,
    modeld, projd, view,
    &winx, &winy, &winz);

  printf("my  %lf %lf %lf\n", winx, winy, winz);
  puts("-----------");

  glMatrixMode(GL_PROJECTION);
  glPopMatrix();

  glMatrixMode(GL_MODELVIEW);
  glPopMatrix();

}

実行結果


glu -269.362487 334.563333 0.978900
my  -269.362487 334.563333 0.978900
-----------
glu -127.939462 262.458218 0.960356
my  -127.939462 262.458218 0.960356
-----------
glu -207.160157 362.861399 0.941500
my  -207.160157 362.861399 0.941500
-----------
glu -63.012933 385.110848 0.985399
my  -63.012933 385.110848 0.985399
-----------
glu -353.644733 582.528452 0.940617
my  -353.644733 582.528452 0.940617
-----------
glu -201.429535 654.202066 0.942501
my  -201.429535 654.202066 0.942501
-----------
glu 72.061917 332.346399 0.977161
my  72.061917 332.346399 0.977161
-----------
glu -10.131515 438.149822 0.985127
my  -10.131515 438.149822 0.985127
-----------

コメントを残す

メールアドレスが公開されることはありません。

日本語が含まれない投稿は無視されますのでご注意ください。(スパム対策)


この記事のトラックバックURL: