2023-05-08(書いているのは2023-11-17)
注!)コードをvisual studioからブログに移す時に間違えてどこかを消してしまったりして正しく動かない場合もあるかもしません。また強度の計算もしていません。
Unityで大気光学現象のシミュレーションはできるのか気になったので挑戦してみる4
「光の気象学※1」という本を買って読んでみると驚くべき事に太陽高度による屈折率の変化の式が載っていたので、これで必要なものは揃いました。
では必要な関数と求める過程を書いていきます。
これは屈折率n、頂角aのプリズムを通る光を表しています。
偏向角Dはiとrの差とi'とr'の差の和なので、
また
rとr'の和はaなので(理由は省きますが簡単です)
これをに代入して
屈折角は、スネルの法則
(n1 は媒質1の屈折率、n2 は媒質2の屈折率)
からrについて解き、n1に空気の屈折率である1に代入して、n2をnにすると
r'は
から
i'はスネルの法則を変形させて、i'について解き、n2に1を代入して、n1をnにすると
これでDを求めることができます。
次に太陽高度による屈折率(n')ですが、光の気象学(柴田清考 著)によると
(hは水平面に対しての角度)
らしいのでそのまま使います。説明を読んでも自分には難しくて分かりませんでした。
臨界角を求める式は、媒質1(屈折率n1)から媒質2(屈折率n2)に光が出ていくとき、
臨界角
今回媒質1の屈折率はnなので、n1をnにし、媒質2は空気なのでn2に空気の屈折率である1を代入して
臨界角
私はこれらをプログラムして解いていくので不要ですが、これらの式をすべてまとめたものがこれです。
(合っているか自信ないですが…)
ではUnityに移ります。
計算を実行するスクリプトmanagerをつくり、空(読み:から)のオブジェクトにアタッチします。
太陽高度の変数h
屈折率の変数n
計算する光線の数の変数number_of_Rays
描画(計算)した数の変数count
これらを宣言し、
using System.Collections; using System.Collections.Generic; using UnityEngine; public class manager : MonoBehaviour {
public float h = 0; public float n = 1.309f; public int number_of_Rays = 100; int count = 0; }
hとnとnumber_of_RaysはUnityから変更できるようにpublicにしておきます。
nには氷の屈折率である1.309※2 を、number_of_Raysにはとりあえず100を代入しておきます。
parhelion(幻日)という名のメソッドをつくってUpdateメソッドで常に呼び出すようにします。
using System.Collections; using System.Collections.Generic; using UnityEngine;
public class manager : MonoBehaviour {
public float h = 0;//太陽高度
public float n = 1.309f;//屈折率
public int number_of_Rays = 100; int count = 0;
void Update() { parhelion(); } void parhelion() { float i = UnityEngine.Random.Range(0.0f, 90.0f);//入射角 float a = 60;//プリズム角 } }
parhelionメソッド内に
入射角の変数i
プリズム角の変数a
を宣言し、iの絶対値は0から90なので、UnityEngine.Random.Range(0.0f, 90,0f)と書きます。
幻日は60°プリズムの現象なのでaには60を代入しておきます。
光が斜めに入った時の屈折率と臨界角を求める部分を書きます。
using System.Collections;
using System.Collections.Generic;
using UnityEngine;
using System; public class manager : MonoBehaviour
{
public float h = 0;//太陽高度
public float n = 1.309f;//屈折率
public int number_of_Rays = 100;
int count = 0;
void Update() { parhelion(); }
void parhelion() { float i = UnityEngine.Random.Range(0.0f, 90.0f);//入射角
float a = 60;//プリズム角
double fh_n = Math.Sqrt(((n * n) - (Math.Sin(h * (Math.PI / 180)) * Math.Sin(h * (Math.PI / 180)))) / (1 - Math.Sin(h * (Math.PI / 180)) * Math.Sin(h * (Math.PI / 180))));//光が斜めに入った時の屈折率
double critical_angle = Math.Asin(1 / fh_n) * 180 / Math.PI;//臨界角 } }
計算でMathクラスを使いたいので
using System;
を追加します。
double b = Math.Asin(a) * 180 / Math.PI;
このように書きます。bの部分の名前は変えても問題ありません。aの部分には数字や変数を入れてください。
arctanやarccosの場合はAsinの部分をそれぞれAtan、Acosと変えればOKです。
もし上手くいかない場合はusing System;を書いてあるか確認してみて下さい。
光が斜めに入った時の屈折率の変数fh_n
臨界角critical_angle
をそれぞれdouble型で宣言して、公式に基づいてコードを書きます。
同じようにrとr'を求めるコードも書きます。ダッシュが変数名として使えなかったので代わりに2を使っています。
using System.Collections; using System.Collections.Generic; using UnityEngine; using System; public class manager : MonoBehaviour
{
public float h = 0;//太陽高度
public float n = 1.309f;//屈折率
public int number_of_Rays = 100;
int count = 0;
void Update()
{
parhelion();
}
void parhelion()
{
float i = UnityEngine.Random.Range(0.0f, 90.0f);//入射角
float a = 60;//プリズム角
double fh_n = Math.Sqrt(((n * n) - (Math.Sin(h * (Math.PI / 180)) * Math.Sin(h * (Math.PI / 180)))) / (1 - Math.Sin(h * (Math.PI / 180)) * Math.Sin(h * (Math.PI / 180))));//光が斜めに入った時の屈折率
double critical_angle = Math.Asin(1 / fh_n) * 180 / Math.PI;//臨界角
double r = Math.Asin(Math.Sin(i * (Math.PI / 180)) / fh_n) * (180 / Math.PI); double r2 = a - r; } }
i'やDも求めていきますが、注意しないといけないのが i の角度によっては r' が臨界角を超えてしまうことがありその時は残りの計算はできなくなるので、ifでr'≦臨界角の場合のみ実行するようにします。
using System.Collections;
using System.Collections.Generic;
using UnityEngine;
using System;
public class manager : MonoBehaviour
{
public float h = 0;//太陽高度
public float n = 1.309f;//屈折率
public int number_of_Rays = 100;
int count = 0;
void Update()
{
parhelion();
}
void parhelion()
{
float i = UnityEngine.Random.Range(0.0f, 90.0f);//入射角
float a = 60;//プリズム角
double fh_n = Math.Sqrt(((n * n) - (Math.Sin(h * (Math.PI / 180)) * Math.Sin(h * (Math.PI / 180)))) / (1 - Math.Sin(h * (Math.PI / 180)) * Math.Sin(h * (Math.PI / 180))));//光が斜めに入った時の屈折率
double critical_angle = Math.Asin(1 / fh_n) * 180 / Math.PI;//臨界角
double r = Math.Asin(Math.Sin(i * (Math.PI / 180)) / fh_n) * (180 / Math.PI);
double r2 = a - r;
if (r2 <= critical_angle) { double i2 = Math.Asin(fh_n * (Math.Sin(r2 * (Math.PI / 180)))) * 180 / Math.PI; double D = i + i2 - a; } } }
i'やDも同じように書いて数式の部分は完成です。
i'はスクリプト内ではi2としています。
試しにDebug.Log(D);を追加して表示される値を見てみると、最小偏角21.7°あたりに集中しており、またそれ以下の数値はないので間違っていないようです。
次はプロットです。描画の方法が分からないので球体のオブジェクトを使います。
オブジェクトを回転させる必要がありますが、UnityのTransformにはdoubleは使えないのでDをfloat型に変換します。
using System.Collections;
using System.Collections.Generic;
using UnityEngine;
using System;
public class manager : MonoBehaviour
{
public float h = 0;//太陽高度
public float n = 1.309f;//屈折率
public int number_of_Rays = 100;
int count = 0;
void Update()
{
parhelion();
}
void parhelion()
{
float i = UnityEngine.Random.Range(0.0f, 90.0f);//入射角
float a = 60;//プリズム角
double fh_n = Math.Sqrt(((n * n) - (Math.Sin(h * (Math.PI / 180)) * Math.Sin(h * (Math.PI / 180)))) / (1 - Math.Sin(h * (Math.PI / 180)) * Math.Sin(h * (Math.PI / 180))));//光が斜めに入った時の屈折率
double critical_angle = Math.Asin(1 / fh_n) * 180 / Math.PI;//臨界角
double r = Math.Asin(Math.Sin(i * (Math.PI / 180)) / fh_n) * (180 / Math.PI);
double r2 = a - r;
if (r2 <= critical_angle)
{
double i2 = Math.Asin(fh_n * (Math.Sin(r2 * (Math.PI / 180)))) * 180 / Math.PI;
double D = i + i2 - a;
float floatD = (float)D; } } }
光のオブジェクトをスクリプトから作りたいのでGameObject型の光のプレハブ用の変数"RayPrefab"を宣言します。
using System.Collections;
using System.Collections.Generic;
using UnityEngine;
using System;
public class manager : MonoBehaviour
{
public float h = 0;//太陽高度
public float n = 1.309f;//屈折率
public int number_of_Rays = 100;
int count = 0;
public GameObject RayPrefab;
void Update()
{
parhelion();
}
void parhelion()
{
float i = UnityEngine.Random.Range(0.0f, 90.0f);//入射角
float a = 60;//プリズム角
double fh_n = Math.Sqrt(((n * n) - (Math.Sin(h * (Math.PI / 180)) * Math.Sin(h * (Math.PI / 180)))) / (1 - Math.Sin(h * (Math.PI / 180)) * Math.Sin(h * (Math.PI / 180))));//光が斜めに入った時の屈折率
double critical_angle = Math.Asin(1 / fh_n) * 180 / Math.PI;//臨界角
double r = Math.Asin(Math.Sin(i * (Math.PI / 180)) / fh_n) * (180 / Math.PI);
double r2 = a - r;
if (r2 <= critical_angle)
{
double i2 = Math.Asin(fh_n * (Math.Sin(r2 * (Math.PI / 180)))) * 180 / Math.PI;
double D = i + i2 - a;
float floatD = (float)D;
}
}
}
Dをそのまま回転の角度に代入すると、iの変域は0~90なのでDは正の数になり、負の領域に描画することができないので、-Dの角度の光も必要です。
まず正の値(右側)を描画するコードを書きます。
using System.Collections;
using System.Collections.Generic;
using UnityEngine;
using System;
public class manager : MonoBehaviour
{
public float h = 0;//太陽高度
public float n = 1.309f;//屈折率
public int number_of_Rays = 100;
int count = 0;
public GameObject RayPrefab;
void Update()
{
parhelion();
}
void parhelion()
{
float i = UnityEngine.Random.Range(0.0f, 90.0f);//入射角
float a = 60;//プリズム角
double fh_n = Math.Sqrt(((n * n) - (Math.Sin(h * (Math.PI / 180)) * Math.Sin(h * (Math.PI / 180)))) / (1 - Math.Sin(h * (Math.PI / 180)) * Math.Sin(h * (Math.PI / 180))));//光が斜めに入った時の屈折率
double critical_angle = Math.Asin(1 / fh_n) * 180 / Math.PI;//臨界角
double r = Math.Asin(Math.Sin(i * (Math.PI / 180)) / fh_n) * (180 / Math.PI);
double r2 = a - r;
if (r2 <= critical_angle)
{
double i2 = Math.Asin(fh_n * (Math.Sin(r2 * (Math.PI / 180)))) * 180 / Math.PI;
double D = i + i2 - a;
float floatD = (float)D;
GameObject Ray_right = Instantiate(RayPrefab);
Ray_right.transform.position = new Vector3(0, 0, 0);
Ray_right.transform.Rotate(-h, floatD, 0); } } }
parhelionメソッド内にプレハブを生成するコード
GameObject Ray_right = Instantiate(RayPrefab);
を書いて
Ray_right.transform.position = new Vector3(0, 0, 0);
Ray_right.transform.Rotate(-h, floatD, 0);
で生成したプレハブを回転させます。-hなのは回転の方向が逆だったからです。
左側は-floatD回転させるので
using System.Collections;
using System.Collections.Generic;
using UnityEngine;
using System;
public class manager : MonoBehaviour
{
public float h = 0;//太陽高度
public float n = 1.309f;//屈折率
public int number_of_Rays = 100;
int count = 0;
public GameObject RayPrefab;
void Update()
{
parhelion();
}
void parhelion()
{
float i = UnityEngine.Random.Range(0.0f, 90.0f);//入射角
float a = 60;//プリズム角
double fh_n = Math.Sqrt(((n * n) - (Math.Sin(h * (Math.PI / 180)) * Math.Sin(h * (Math.PI / 180)))) / (1 - Math.Sin(h * (Math.PI / 180)) * Math.Sin(h * (Math.PI / 180))));//光が斜めに入った時の屈折率
double critical_angle = Math.Asin(1 / fh_n) * 180 / Math.PI;//臨界角
double r = Math.Asin(Math.Sin(i * (Math.PI / 180)) / fh_n) * (180 / Math.PI);
double r2 = a - r;
if (r2 <= critical_angle)
{
double i2 = Math.Asin(fh_n * (Math.Sin(r2 * (Math.PI / 180)))) * 180 / Math.PI;
double D = i + i2 - a;
float floatD = (float)D;
GameObject Ray_right = Instantiate(RayPrefab);
Ray_right.transform.position = new Vector3(0, 0, 0);
Ray_right.transform.Rotate(-h, floatD, 0);
GameObject Ray_left = Instantiate(RayPrefab);
Ray_left.transform.position = new Vector3(0, 0, 0);
Ray_left.transform.Rotate(-h, -floatD, 0);
} } }
右と左で2個のドットを描画したのでcountを2増やします。
using System.Collections;
using System.Collections.Generic;
using UnityEngine;
using System;
public class manager : MonoBehaviour
{
public float h = 0;//太陽高度
public float n = 1.309f;//屈折率
public int number_of_Rays = 100;
int count = 0;
public GameObject RayPrefab;
void Update()
{
parhelion();
}
void parhelion()
{
float i = UnityEngine.Random.Range(0.0f, 90.0f);//入射角
float a = 60;//プリズム角
double fh_n = Math.Sqrt(((n * n) - (Math.Sin(h * (Math.PI / 180)) * Math.Sin(h * (Math.PI / 180)))) / (1 - Math.Sin(h * (Math.PI / 180)) * Math.Sin(h * (Math.PI / 180))));//光が斜めに入った時の屈折率
double critical_angle = Math.Asin(1 / fh_n) * 180 / Math.PI;//臨界角
double r = Math.Asin(Math.Sin(i * (Math.PI / 180)) / fh_n) * (180 / Math.PI);
double r2 = a - r;
if (r2 <= critical_angle)
{
double i2 = Math.Asin(fh_n * (Math.Sin(r2 * (Math.PI / 180)))) * 180 / Math.PI;
double D = i + i2 - a;
float floatD = (float)D;
GameObject Ray_right = Instantiate(RayPrefab);
Ray_right.transform.position = new Vector3(0, 0, 0);
Ray_right.transform.Rotate(-h, floatD, 0);
GameObject Ray_left = Instantiate(RayPrefab);
Ray_left.transform.position = new Vector3(0, 0, 0);
Ray_left.transform.Rotate(-h, -floatD, 0);
count += 2; } } }
あとは太陽も描画したいので、太陽用のGameObject型の変数SunPrefabをつくります。
Startメソッドを追加して、最初に描画するようにします。
using System.Collections;
using System.Collections.Generic;
using UnityEngine;
using System;
public class manager : MonoBehaviour
{
public float h = 0;//太陽高度
public float n = 1.309f;//屈折率
public int number_of_Rays = 100;
int count = 0;
public GameObject RayPrefab;
public GameObject SunPrefab; void Start() { GameObject Sun = Instantiate(SunPrefab); Sun.transform.position = new Vector3(0, 0, 0); Sun.transform.Rotate(-h, 0, 0); } void Update()
{
parhelion();
}
void parhelion()
{
float i = UnityEngine.Random.Range(0.0f, 90.0f);//入射角
float a = 60;//プリズム角
double fh_n = Math.Sqrt(((n * n) - (Math.Sin(h * (Math.PI / 180)) * Math.Sin(h * (Math.PI / 180)))) / (1 - Math.Sin(h * (Math.PI / 180)) * Math.Sin(h * (Math.PI / 180))));//光が斜めに入った時の屈折率
double critical_angle = Math.Asin(1 / fh_n) * 180 / Math.PI;//臨界角
double r = Math.Asin(Math.Sin(i * (Math.PI / 180)) / fh_n) * (180 / Math.PI);
double r2 = a - r;
if (r2 <= critical_angle)
{
double i2 = Math.Asin(fh_n * (Math.Sin(r2 * (Math.PI / 180)))) * 180 / Math.PI;
double D = i + i2 - a;
float floatD = (float)D;
GameObject Ray_right = Instantiate(RayPrefab);
Ray_right.transform.position = new Vector3(0, 0, 0);
Ray_right.transform.Rotate(-h, floatD, 0);
GameObject Ray_left = Instantiate(LayPrefab);
Ray_left.transform.position = new Vector3(0, 0, 0);
Ray_left.transform.Rotate(-h, -floatD, 0);
count += 2;
} } }
最後にCount≦number of raysの時のみparhelionメソッドを呼び出すようにします。
using System.Collections;
using System.Collections.Generic;
using UnityEngine;
using System;
public class manager : MonoBehaviour
{
public float h = 0;//太陽高度
public float n = 1.309f;//屈折率
public int number_of_Rays = 100;
int count = 0;
public GameObject RayPrefab;
public GameObject SunPrefab;
void Start()
{
GameObject Sun = Instantiate(SunPrefab); Sun.transform.position = new Vector3(0, 0, 0);
Sun.transform.Rotate(-h, 0, 0);
}
void Update()
{
if (count < number_of_Rays) { parhelion(); } }
void parhelion()
{
float i = UnityEngine.Random.Range(0.0f, 90.0f);//入射角
float a = 60;//プリズム角
double fh_n = Math.Sqrt(((n * n) - (Math.Sin(h * (Math.PI / 180)) * Math.Sin(h * (Math.PI / 180)))) / (1 - Math.Sin(h * (Math.PI / 180)) * Math.Sin(h * (Math.PI / 180))));//光が斜めに入った時の屈折率
double critical_angle = Math.Asin(1 / fh_n) * 180 / Math.PI;//臨界角
double r = Math.Asin(Math.Sin(i * (Math.PI / 180)) / fh_n) * (180 / Math.PI);
double r2 = a - r;
if (r2 <= critical_angle)
{
double i2 = Math.Asin(fh_n * (Math.Sin(r2 * (Math.PI / 180)))) * 180 / Math.PI;
double D = i + i2 - a;
float floatD = (float)D;
GameObject Ray_right = Instantiate(RayPrefab);
Ray_right.transform.position = new Vector3(0, 0, 0);
Ray_right.transform.Rotate(-h, floatD, 0);
GameObject Ray_left = Instantiate(RayPrefab);
Ray_left.transform.position = new Vector3(0, 0, 0);
Ray_left.transform.Rotate(-h, -floatD, 0);
count += 2;
}
}
}
これで計算側のスクリプトは完成です。
Unityに戻って光と太陽のオブジェクトをつくります。サイズは精密な計算をしたわけではなく適当に決めました。
これらをプレハブにし、そのプレハブを空のオブジェクトにアタッチしたmanagerスクリプトのRayPrefabとSunPrefabの欄にドロップします。
今の状態で実行しても、すべての光と太陽のプレハブは(0,0,0)に集まっているので、一定の距離進ませないといけません。
一定の距離でとめるにはどうすれば良いでしょうか。for文でするのも考えましたが後で色を変更するためややこしくなりそうです。なので原点中心の巨大な球体を置き、それとの衝突が終わったら止めるという方法にします。
これは光のオブジェクトの方に書いていきます。
光を管理する用のスクリプトRayをつくり、書いていきます。
これはSphereタグの付いたトリガーのオブジェクトとの衝突が終わったら…というコードです。
using System.Collections;
using System.Collections.Generic;
using UnityEngine;
public class Ray : MonoBehaviour
{
void OnTriggerExit(Collider other) { if(other.gameObject.tag == "Sphere") { } } }
bool型のIsStopという変数をつくり、falseを代入します。
OnTriggerExitメソッド内のif文の中にIsStop = true;を追加し、衝突が終わったらIsStopがtrueになるようにします。
using System.Collections;
using System.Collections.Generic;
using UnityEngine;
public class Ray : MonoBehaviour
{
bool IsStop = false;
void OnTriggerExit(Collider other)
{
if(other.gameObject.tag == "Sphere")
{ IsStop = true; } } }
Updateメソッドを追加し、if文でIsStopがfalseの時のみz方向に1ずつ動かすようにします。
using System.Collections;
using System.Collections.Generic;
using UnityEngine;
public class Ray : MonoBehaviour
{
bool IsStop = false;
void Update()
{
if(IsStop == false) { this.transform.Translate(0, 0, 1.0f); } }
void OnTriggerExit(Collider other)
{
if(other.gameObject.tag == "Sphere")
{
IsStop = true;
}
} }
これで一定の距離で止める事ができるようになりました。
Unityに戻って巨大な球体のオブジェクトを置き、Sphere ColliderのIs Triggerにチェックを入れます。OnTriggerExitを使うので、Rigidbodyコンポーネントを追加してUse Gravityのチェックを外し、重力が作用しないようにします。そしてSphereという名前のタグをつくって付けます。
一色の背景の方が見やすいので、Blenderで中が空洞の球体をつくり、それをエクスポートしてUnityに持ってきます。さっき置いた球体以上の大きさで、尚且つ描画される範囲に収まるならどんな大きさでもよいです。光を考慮したくないので、マテリアルをつくってShaderをUnlit/Colorにしてとりあえず白色にします。
背景を白色にしたので、光と太陽の色を白から黒色に変た方が良いのですが、変えると生成されて止まるまでの間、飛んでいく球体が見えてしまいます。(下の動画参考。#1の時の様子です) 邪魔なので、止まった時に表示するようにしようと思います。
RayスクリプトにMaterial型の変数を二つ宣言して、それぞれfirstColorとSecondColorという名前にします。
firstの方は動いているとき用の色で、secondは止まった時用の色です。それぞれIsStopがfalseの時とそうでない時の所にその色に変更するコードを書きます。
using System.Collections;
using System.Collections.Generic;
using UnityEngine;
public class Ray : MonoBehaviour
{
bool IsStop = false; public Material firstColor; public Material SecondColor;
void Update() {
if(IsStop == false)
{
this.transform.Translate(0, 0, 1.0f); GetComponent<Renderer>().material.color = firstColor.color; } else { GetComponent<Renderer>().material.color = SecondColor.color; } }
void OnTriggerExit(Collider other)
{
if(other.gameObject.tag == "Sphere")
{
IsStop = true;
}
}
}
Unityに戻り、マテリアルをつくります。firstのマテリアルは見えないように透明にしておきます。透明度を変更するためにRendering ModeをFadeにしておきます。
secondの方は適当に透明度が約16%(=40/255)ほどの黒色にしておきます。
太陽用のsecondも作ります。こちらは半透明である必要はないので、Shaderを反射を防げるUnlit/Colorを使い、真っ黒にします。
光のプレハブのRayスクリプトにはこのように↓マテリアルを入れます
太陽は↓(Sunは太陽用の不透明のマテリアル)
これで完成です。
太陽高度0°で実行すると、うまく幻日が表示されました。
強度の計算をしていないので端の方は実際より強く表示されているはずです。
太陽高度が上がるにつれて幻日が離れていくのも再現されています。下から二つ目では太陽に近くなっていますが、レンズの影響です。
幻日は太陽高度が61°弱で消えますが、このシミュレーションでも60.9°で消えています。
光線の数を6000ほどにしても重くならなかったので、これからこの方式でしていくことになりそうです。
式を書き始めてから描画できるまで4日かかりました。
今はまだ揺れや強度などには対応していないのでこれからそれらも取り入れていきたいです。