旋回する。（フライトシミュ）

[Ma]

＞揚力のベクトルを求める。（フライトシミュ）
http://www.play21.jp/board/formz.cgi?ac ... &rln=48360
のスレでお世話になりました。
おかげさまで無地に実装できました。


揚力実装後に新たな点でつまづいたので質問させていただきます。
よし、これできれいに旋回できるかなーなんておもってたんですが、そんなわけありませんでした。

実際の飛行機は操縦桿の操作（ひいたり、おしたり）し、（水平尾翼の）揚力を増減によって機体の状態角度を変ます。
この点は、今の私のゲームではキーボードの入力から直接角度を変えているため、同様な操作ができます。

この後が問題です。
実世界での飛行機は、機体のこの状態角度の変化によって風がより主翼にあたり、揚力がさらに上昇し旋回を手助けます。
ですが、私のゲームでの現状は、状態角度を変えただけで揚力の増減が発生しません。
結果的に、移動している方向と機体が向いている方向がぜんぜんかけはなれてしまいました。

では、揚力に係数をかけてあげればよくて、その係数は機体の面と風が吹く方向によって決まるはずだ。

ということで、
まず、機体の進行方向のベクトルと風のベクトルを外積してみましたが、この方法だとロール角が考慮されないので、これだけでは駄目でした。
（ロール角によっては風と主翼の面がまったくあたらないような角度もあるため。）

じゃぁ、外積したあとに、さらにロール角を使って何か計算するだろう。
と思うまではたどり着きましたが、どうロール角を扱えばいいのか分かりません。

これが質問です。



ちょっと違うアスペクトから言いますと、求めたいのはある位置（機体のベロシティの単位位置ベクトル）から、原点に置かれ任意の角度に回転させた正方形を見たときの、その面（機体の主翼）の面積のような値です。
ただし、裏面が見えたときはプラスで、表面が見えたときはマイナスになるような値。

『このような特徴』を持った値は、揚力の係数として扱うことができると考えています。
（正方形の話は、あくまでイメージを伝えるために書きました。）

よろしくおねがいします。

[softya]

だんだん、真面目なフライトシム化してますね。

NekoFlightの事を思い出したので、参考に成るかと貼っておきます。
http://homepage1.nifty.com/kaneko/nflight2.htm
http://homepage1.nifty.com/kaneko/nflightm.htm
ちなみに、Winnyで有名な金子勇氏のページです。
色々ここには面白いものが置いてあります。
http://homepage1.nifty.com/kaneko/index.htm

とりあえず「進行方向のベクトルと風のベクトルを外積」は方向性として違うと思いますよ。

[Ma]

「翼計算」
のところを読んでいたら参考になるようなことがありました。
いつも的確で短いご指摘、感謝してます。

さっそく、手始めとして迎角を求めようと思います。
面と線との間の角ですね。

ちょっと調べてみます。

追記
先ほど求めた（翼に対する法線ベクトル）揚力ベクトルから、ベロシティのベクトルとの内積からその角度をもとめる。
９０度からそれをひいて、面と線の角度ですね。
表面か裏面かもわかりますね。

あとは例の係数と方程式にのっとって・・・と。
いったん休憩にしますｗ


さらにPS
そういえば、DXライブラリ3D のほうで、フライトシミュレーター作っていた方がいたとおもったんだけど、どうなったんだろう。。。。

[Ma]

なんとか応用できましたが、設計からいろいろひどい部分がでてきたので、
しばらく修正していそうです。
とりあえず、質問はいただいたページにかかれていることを実装してからにします。

[Ma]

どうも、うまくいきませんでした。。。orz
まだまだ、力のモーメントを考える前の段階です。

旋回はなんとなくできるようになってきましたが、
やはりまだまだ異常な動きがあります。


特になんとかしたいのは、機首をy軸方向にゆっくりむけながら上昇すると
速度がいきなり暴走し始めます。

なんだか、理由がよくわかりそうにないので、よければソースを読んでいただけるとうれしいです。

[softya]

とりあえず、じっくり見たいので夜まで待って下さいね。
起動しみたら動かず(汗)

[Ma]

お手数おかけします＞＜


http://www1.axfc.net/uploader/Sc/so/85498
PASS: sim

↑　いまUPした、ちょい改良版もあるので、どちらを見ていただいてもかまいません。
どちらとも同じような問題をかかえてます。。。

速度暴走というかなんというか、、、ともかく暴走しますｗｗｗ

キー操作：
↑キー　すらすとあげ
↓キー　すらすとさげ
Wキー　えれべーたーあげ（ある程度速度でないと意味なし）
Sキー　えれべーたーさげ
Dキー　みぎロール（地上では効果なし）
Aキー　ひだりろーる（地上では効果なし）

QとEを　地上移動方向＆らだー　をつけようと思いましたが、まだできてません。
最初に滑走しはじめて機首あげればとびます。

[softya]

すいません、ちょっと時間が取れないので明日にさせて下さい。
出来れば簡単に今やっている処理を教えて貰うと解読時に助かります。

[Ma]

とりあえず、今行っているアルゴリズムは以下のような流れです。
void player::update()
内での話しになります。

以下、空気密度、とか翼面積だとかでてきますが、プレイ中は定数です。

１．翼対して、垂直な単位法線ベクトルを得る。（機体の進行方向とは関係なく、機首の向きに垂直）（クォータニオンを使う）
２. 揚力方向の単位ベクトルを求める。（法線とは違い、進行方向に対して垂直なベクトル。）（クォータニオンを使う）
３．getMukaeKaku 関数に法線を渡し、翼と進行方向の間の角度を求める。
４．公式　揚力＝揚力係数＊翼速度＾2＊(翼面積＊空気密度／2)　　と、 揚力係数 ∝ 3で得た迎角　を利用して、2　で得た揚力ベクトルの長さを求める。
５．揚力ベクトルを速度ベクトルに加算。
６．公式　抗力＝抗力係数＊翼速度＾2＊(翼面積＊空気密度／2)　と　抗力係数 ∝ mukaeKaku^2 を利用して、抗力ベクトルの長さを求める。
７．速度ベクトルに抗力ベクトルを加算。
８．重力を加算。
９．推力を加算。
１０．エレベーターの角度elevator と（クォータニオンを使う）getRotatedVec関数を使って、エレベーターの迎角を求める。
１１．エレベーターによって、機首角度がかわる事を模擬的に計算する。（力のモーメントが面倒だったので、今は模擬的になっている）
（ちなみに、エレベーターによる揚力の角速度(模擬的)∝揚力係数＊翼速度＾2＊(翼面積＊空気密度／2)）
１２. エレベーターによる揚力があるとき、ロール角があったら仰角だけでなく旋回することもあるので、rotAng によって成分分けする。
１３．UI処理と、角度調節。
１４．速度を位置に加算。
（１５．カメラの位置、角度設定と情報表示。）

図解を添付します。

[Ma]

追記

垂直尾翼に対して、まだ計算していないので
すらすとを切って、旋回して水平に戻すと機体がスライド飛行しますｗ


離陸してしばらくは、とても良好な動きをしてくれるんですが、しばらくすると、とんでもないことになるんですよね。。。

[softya]

軽く見ましたが、フライトシムを真面目に作った事が無いので勘が働きませんね(^^;
っことで提案なのですが、飛ばす前に風洞実験みたいに翼の角度と風で揚力ベクトルがどう働いているか実験施設のシミュレータみたいな物を作ってみてはどうでしょう？機体は動かずに一点に固定で、機体の回転や風速だけ変えれるって本当に風洞実験の施設みたいにするわけですね。
各ベクトル値は、視覚で分かるように矢印か棒をベクトル毎に色分けしてポリゴンで表示してやると異常値が出たら視覚的にすぐ分かると思います。
あとは、ラダーやエレベータも風洞実験して搭載していくって形がわかり易いと思うんですが。

>離陸してしばらくは、とても良好な動きをしてくれるんですが、しばらくすると、とんでもないことになるんで
すよね。。。
この問題も風洞実験で分かるかと思いますよ。

[Ma]

おぉ、それはいいアイディアですね！
次に時間ができたときさっそくやってみますね。

[Ma]

空洞実験のようなものをやっていろいろ問題がみつかり、大半は解決しましたが、、、
新しい問題にぶちあたりました。

どうやら、進行方向が縦に９０度を超える瞬間に揚力の向きが反転してしまうことがわかりました。
原因は簡単で、移動方向のｘｚ座標上の角度を
float hMoveAngle = atan2(myVel.x,myVel.z);
と、していたためです。
９０度に越えた瞬間XZ座標上の水平角度（hMoveAngle）が180度反転したため、機体がy軸を対象に１８０度回転したようになり、当然揚力も反転します。
結果、９０度に達した瞬間に（水平角度１８０度反転を）いったりきたり毎フレームごとに繰り返すことになりました。


法線ベクトルの時は、機首の方向（UIによって変化する）から導きだすので、鉛直角度の範囲は－∞度から∞度でした。
（実際は-180度から180度に制限しますが。）

が、揚力ベクトルは進行方向（機体の状態によって様々に変化）とatan2関数から導き出すので、鉛直範囲がー９０度から９０度になってしまう。

なんとか -180度から180度までわかるようにするには、前フレームの鉛直角度から今フレームの鉛直角度の変化がわかればいいとおもったんですが、
（８９度から９１度にはいる瞬間がわかればいいとおもったから）
その８９度の次の瞬間に、機首をおろしたのか、それともあげつづけたのかは分からないので、
９２度と８７度の違いが区別できないことになります。
だからといって、機首の操作だけで場合わけするわけにはいきません。
ラダーだとか他のさまざまな操作を使った結果に反転したのかもしれないためです。

この問題、なんとかならないでしょうか？

質問をまとめると、この進行方向の仰角の範囲を -90度から90度を、-180度から180度までになおすことはできませんか？
よろしくおねがいします。

[Ma]

↓試してみたこと。
この方法でなんとかいけるか、とおもったんですが、どうも正しく動作しません。。。。
プロジェクトも添付します。
http://www1.axfc.net/uploader/Sc/so/86360
PASS　　sim

VECTOR getLift(){

float size = sqrt(myVel.z*myVel.z+myVel.x*myVel.x);

//printfDx("%4.2f\n",hMoveAngle);
float TempVerticalMoveAngle = -atan2(myVel.y,size);
float hMoveAngle = atan2(myVel.x,myVel.z);

static float vMoveAngle = TempVerticalMoveAngle;
static float PreHMoveAngle = hMoveAngle;
static float PreVMoveAngle = TempVerticalMoveAngle;




float vAngleDifferece = vMoveAngle-PreVMoveAngle;
float hAngleDifferece = hMoveAngle-PreHMoveAngle;
PreVMoveAngle = vMoveAngle;
PreHMoveAngle = hMoveAngle;

static bool vAngInverted = false;
if(abs(abs(hAngleDifferece)-PHI_F) < 0.04f){//反転
vAngInverted = !vAngInverted;
printfDx("flipped:%4.2f %4.2f\n",vMoveAngle);
}



if(vAngInverted){
// printfDx("vMoveAngle%4.2f,ho%4.2f diff%4.2f diff%4.2f",vMoveAngle,hMoveAngle,vAngleDifferece,hAngleDifferece);
vMoveAngle += vAngleDifferece;
hMoveAngle += PHI_F;
// printfDx("vMoveAngle%4.2f,ho%4.2f\n",vMoveAngle,hMoveAngle);

}
else{
vMoveAngle = TempVerticalMoveAngle;
}



DrawFormatString(0,40,GetColor(255,255,255),"MovingAng(V%4.2frad,H%4.2frad)",abs(vMoveAngle)<0.01f?0.0f:vMoveAngle,hMoveAngle);
return getRotatedVec(hMoveAngle,vMoveAngle,rotAng,PHI_F/2);
}




キー追加操作
I カメラあげ
K カメラさげ
J カメラひだり
L かめらみぎ
U カメラずーむいん
O　かめらずーむあうと

[softya]

またまた、ざっと見ただけで申し訳ないのですが、進行方向ベクトルに対して垂直なベクトルがｙ方向にプラスかマイナスかで背面飛行と正面？飛行が区別できると思いますので、それでatan2の結果を-180から+180に補正できませんでしょうか？

残念ながら私のところでちゃんと動かないんですよね。私のノートＰＣの問題だと思うんですが。

[Ma]

＞進行方向ベクトルに対して垂直なベクトルがｙ方向にプラスかマイナスかで背面飛行と正面？飛行が区別できると

進行方向に垂直なベクトル（かつ機体↑方向）なのは、「揚力」で、この揚力を求めるための計算（クォータニオン）をするために進行方向のXZ平面との仰角を求めようとしています。

つまり、
進行方向とXZ平面との仰角→揚力ベクトル
しようとしているのに
揚力ベクトル→進行方向とXZ平面との仰角　はできません。
なので、揚力（進行方向ベクトルに対して垂直なベクトル）を使うのはできないんじゃないかな、とおもいます。

くりかえしになりますが、進行方向ベクトルと機首の向きベクトルは異なります。
進行方向は単純に揚力などの力によって得られたベクトル（X,Y,Zの値）ですが、機首のベクトルはユーザーインターフェイスから得た角度（XZ平面角度,XZ平面迎角度(-PI～PI)）から算出しています。
この違いが９０度以上の角度を得られない原因になっているんですよね。。。


質問をくりかえすかたちになりますが、進行方向のベクトルは atan2(y方向の速さ,(x,z平面上の速さ)) から求めているので９０度以上になるとUターンしてもどってきます。（８９度→９１度となっても、計測上は８９度→８９度で、平面上の角度は反転する。）

これを、８９度→９１度にする方法が、なかなかおもいつきません。


使える値
・機首の向きベクトルと、角度（迎角については-PI～PI）。（進行方向とは無関係）
・機体の状態の法線ベクトル（進行方向とは無関係）
・進行方向の角度の変化量（ただし、９０度付近はあやふやになる。）
など。


計算しようとしている値
・揚力ベクトル

それに必要な値
・クォータニオンをするために必要な進行方向とXZ平面との仰角（-PI～PI）。


質問
「それに必要な値」の仰角が９０度以上を示さないため、実際の角度が９０度を超えると揚力が反転する。

[softya]

あれ？
進行方向と機体方向は別じゃ無かったでしたっけ？
揚力は、機体方向で決まるのでは？
私が混乱してるのかな？

それはそれとして、垂直なベクトルはちょっとした工夫で外積で求められると思います。
進行方向ベクトルを(x,y,z)とした時に、xz平面に投影したベクトルをまず求めます。(x,0,z)ですね。この(x,0,z)を更に９０度時計方向に回します。９０度回すと(-z,0,x)となります。これをベクトル２とします。
進行方向ベクトルとベクトル２を外積すると、垂直なベクトルを求めることが出来ます。

[Ma]

追記
>進行方向と機体方向は別じゃ無かったでしたっけ？
>揚力は、機体方向で決まるのでは？

>進行方向ベクトルと機首の向きベクトルは異なります。
この違いが揚力を発生させます。
揚力は、風の向きと翼の向きの違いで発生するものだと考えています。
つまり、移動方向と機首方向の間の角度で変動します。



＞それはそれとして、垂直なベクトルはちょっとした工夫で外積で求められると思います。
たしかにそのとおりですが、それは（ロール角によって）平面がねじれてない場合のみです。
それでは、ロール角が考慮されなくなってしまいます。
だからこそのクォータニオンでした。


>計算しようとしている値
>・揚力ベクトル

揚力ベクトルはロール角によって著しく変化します。



ちなみに、さきほど、ようやく１ステップ進みまして、
あまり使いたくなかった方法でしたが、以下の方法で９０度以上の判定がなんとかなりました。
ただ、やり方的に誤作動がないとは言えないので、もっといい方法おもいついたらよろしくおねがいします。


//揚力（進行方向(機首方向ではない)に対して垂直上向きのベクトル。）
VECTOR getLift(){

//ベロシティXZ平面成分
float size = sqrt(myVel.z*myVel.z+myVel.x*myVel.x);

//迎角
float TempVerticalMoveAngle = -atan2(myVel.y,size);
//XZ平面水平角度
float hMoveAngle = atan2(myVel.x,myVel.z);

//垂直角度
float vMoveAngle = TempVerticalMoveAngle;

//前の角度
static float PreHMoveAngle = hMoveAngle;


//前フレームとの水平角度の違い
float hAngleDifferece = hMoveAngle-PreHMoveAngle;
//水平角度記録
PreHMoveAngle = hMoveAngle;


//±９０度を超えたかどうか記録。
static bool flip = false;

//水平角度がかなりかわった→反転した可能性大→揚力を反転する。
if(abs(hAngleDifferece) > PHI_F/2){
flip = !flip;
printfDx("flip %d %d %d\n",(int)(PreHMoveAngle*180/PHI_F),(int)(hMoveAngle*180/PHI_F),(int)(hAngleDifferece*180/PHI_F));
}

//反転しているなら
if(flip){
//迎角修正。
vMoveAngle = -PHI_F-TempVerticalMoveAngle;
//水平反転
hMoveAngle -= PHI_F;
}

return getRotatedVec(hMoveAngle,vMoveAngle,rotAng,PHI_F/2);
}





ただ、新しい問題がまた発生しました＾＾；
というのは、getRotatedVec関数はクォータニオンを利用する関数ですが、abs(vMoveAngle) > PHI_F/2 の時に対応しておらず、おかしな挙動をすることがわかりました。
先にできることをやってから報告したかったのですが、1,2日ほど忙しくなるので先に報告だけさせていただきました。

[softya]

どうやら勘違いしていたようなので説明とプログラムを見直しますね。
ただ、すぐ時間が取れないので、お待ちください。

調べる前に、ひとつだけ聞きたいのは揚力は進行方向に垂直なベクトルと定義されていますが、垂直ベクトルはどの平面に対して垂直なのでしょうか？それが定義出来れば、垂直ベクトルを求められると思うのですが。

[Ma]

お返事おくれました、すいません。
>垂直ベクトルはどの平面に対して垂直なのでしょうか？それが定義出来れば、垂直ベクトルを求められると思うのですが。

垂直ベクトル＝揚力の単位ベクトル　と解釈しますが、
進行（移動）方向のベクトルと期待のロール角からできる面に対して垂直なのが垂直ベクトルです。

それにたいして、ここでいう法線ベクトルは、翼の面に対して垂直なベクトルです。


このへんは、プログラムが動けば一番わかりやすいとおもうのですが、どうしても動かないのでしょうか？
よろしければ、リリースモードでコンパイルして起動できる状態のものをご用意させていただこうかと思います。

[Ma]

誤字
期待→機体

[Ma]

最新版です。
よろしければこちらをごらんになってください。
また、起動できるようにexeファイルも入っています。

PW sim
http://www1.axfc.net/uploader/Sc/so/87342

[softya]

exe版起動出来ました。確かに、妙な機動をしますね。
atan2辺りが問題なのも視覚的に納得出来ました。
ただ、時間が無いのでごめんなさい今日は無理かも知れません。

[Ma]

softya さん、いそがしいところすいません＞＜
それからありがとうございます。



ところで、さきほど原因が発覚しました＾＾；
↓誤り
//回転軸を機体の進行方向に指定し、ロール角だけ回転する時に利用するクォータニオンを生成。
qqq = MakeRotationalQuaternion(-rotAngle, sin(hAngle), sin(-vAngle), cos(hAngle));
rrr = MakeRotationalQuaternion(rotAngle, sin(hAngle), sin(-vAngle), cos(hAngle));




↓正しい
//回転軸を機体の進行方向に指定し、ロール角だけ回転する時に利用するクォータニオンを生成。
qqq = MakeRotationalQuaternion(-rotAngle, sin(hAngle)*cos(vAngle), sin(-vAngle), cos(hAngle)*cos(vAngle));
rrr = MakeRotationalQuaternion(rotAngle, sin(hAngle)*cos(vAngle), sin(-vAngle), cos(hAngle)*cos(vAngle));




進行方向にしないといけないのに、なっていなかったみたいです。
クォータニオン作成が間違っていたため、反転後におかしいことになったみたいです。
たまたま反転する前にロール角で機体を上下反転させたら法線がおかしなことになっていたところから、
きづくことができました。

ここが一番致命的になっていまして、もう一箇所も改良したところ、綺麗にとびました！！
ありがとうございました。
あとは、垂直尾翼の計算をするだけですね。

softya さんのおかげで、モチベーションをキープすることができました。ありがとうございます＞＜

[softya]

うまく動いて良かったです。
あんまりお手伝い出来ませんでしたが( ^ ^ ;
出来上がったら、遊ばせて下さいね。
・・・未だに自分の環境でコンパイルしたものが動かない原因は謎・・・。