<!--
// 天文計算関係
d2r = Math.PI / 180.0;

// 真行動傾斜角の計算（答えはラジアン）
function GetEpci(T)
{
	var Epci;
	T /= 36525.0;
	Epci = 0.00256 * Math.cos((1934 * T + 235) * d2r) + 0.00015 * Math.cos((72002 * T + 201) * d2r);
	Epci += 23.43928 - 0.01301 * T;
	return Epci * d2r;
}

// グリニッジ視恒星時の計算答えは、Hour
function GetGMT(UT,T)
{
	var gmt;
	gmt = UT + 6.69736 + 24.000513 * T / 365.25;
	gmt += 0.00029 * Math.sin((1934 * T / 36525 + 235) * d2r);
	gmt = gmt / 24.0;
	gmt = (gmt -  Math.floor(gmt)) * 24.0;
	return gmt;
}

// 極座標（距離、経度、緯度）から直交座標への変換（角度はラジアン）
function POL2XYZ(p,x)
{
	with(Math)
	{
		x[2] = p[0] * sin(p[2]);
		x[1] = p[0] * cos(p[2]);
		x[0] = x[1] * cos(p[1]);
		x[1] *= sin(p[1]);
	}
}

// 直交座標から極座標への変換
function XYZ2POL(x,p)
{
	var r;
	r = x[0] * x[0] + x[1] * x[1];
	with(Math)
	{
		p[0] = sqrt(x[2] * x[2] + r);
		r = sqrt(r);
	
		p[2] = atan2(x[2],r);
		p[1] = atan2(x[1],x[0]);
	}
}

// 直交座標の回転
// x[0,1,2]をbase軸（0 = x,1 = y, z = 2) の周りに p 回転させる。結果は xに戻す
function ROTXYZ(x,p,base)
{
	var nx = new Array(3);
	var e1,e2;

	if (base == 0) {e1 = 1; e2 = 2;}
	else if (base == 1) {e1 = 0; e2 = 2;}
	else {e1 = 0; e2 = 1;}

	with(Math)
	{
		nx[base] = x[base];
		nx[e1] = cos(p) * x[e1] - sin(p) * x[e2];
		nx[e2] = sin(p) * x[e1] + cos(p) * x[e2];
	}
	for (e1 in nx) x[e1] = nx[e1];
}

// 黄道極座標を赤道極座標に変換する。FLAG == 1 の時黄道->赤道変換。それ以外は、逆変換
// Kou に元の座標値。Sekに戻り値が入る
function Koudou2Sekidou(Kou,Sek,T,FLAG)
{
	var e;
	var xyz = new Array(3);
	e = GetEpci(T);	// 真黄道傾斜角
	if (FLAG != 1) e = -e;	// FLAG != 1 ならば、赤道->黄道変換
	
	POL2XYZ(Kou,xyz); // 黄道座標 -> 赤道座標変換
	ROTXYZ(xyz,e,0);
	XYZ2POL(xyz,Sek);
}


// 地平高度角計算 入力、出力とも角度の単位はラジアンとする。
var Now_dcd;	// 方位角計算に流用するため大域変数にする
function GetEL(T,ra,dc,lng,lat)
{
	var hangl,el,latd,hour;

	hour = T - 0.5;
	hour = (hour - Math.floor(hour)) * 24.0;

	hangl = GetGMT(hour,T) * 15.0 * d2r; // グリニッジ視恒星時計算(ラジアン)
	hangl = ra - (hangl + lng); // 時角
	Now_dcd = 90 * d2r - dc;	// 高度角計算
	latd = 90 * d2r - lat;
	with (Math)
	{
		el = cos(latd) * cos(Now_dcd) + sin(latd) * sin(Now_dcd) * cos(hangl);
		el = atan2(sqrt(1.0 - el * el),el);
	}
	return Math.PI/2.0 - el;
}

// 方位、高度を得る。結果は azel[] にaz,elの順に入る。内部で GetEL() を呼び出す
// 入力、出力とも角度はすべてラジアン単位とする
function GetAZEL(T,ra,dc,lng,lat,azel)
{
	var eld;

	azel[1] = GetEL(T,ra,dc,lng,lat);
	eld = Math.PI/2.0 - azel[1];

	with(Math)
	{
		azel[0] = (cos(Now_dcd) - cos(90 * d2r - lat) * cos(eld));
		azel[0] /= sin(90 * d2r - lat) * sin(eld);
		azel[0] = atan2(sqrt(1.0 - azel[0] * azel[0]),azel[0]);
	}
}

// 方位角計算。内部では、GetAZEL() を呼び出す。
function GetAZ(T,ra,dc,lng,lat)
{
	var azel = new Array(2);
	
	GetAZEL(T,ra,dc,lng,lat,azel);
	return azel[0];
}

// 大気差計算 視高度を引数として補正する大気差を返す。答えは(度)
// -1.0,-0.5は天体位置表の値。あとは理科年表の値
function GetRadau(ah) {
	var y = new Array(56.4583,45.0,34.4,28.65,24.2833,18.2167,14.3333,11.6833,9.8167,7.3667,5.2833,2.6333,1.15,0.55,0);
	var x = new Array(-1.0,-0.5,0.0,0.5,1.0,2.0,3.0,4.0,5.0,7.0,10.0,20.0,40.0,60.0,90.0);
	var h,th,i,j,ns,ne;

	for (ne = 0; ne < y.length ; ne++) {
		if (x[ne] >= ah) break;
	}
	if (ne >= y.length) {
		ne = y.length - 1;
	}
	if (ne < 3) ne = 3;
	ns = ne - 3;
	h = 0.0;
	for (i = ns; i <= ne ; i++) {
		th = y[i];
		for (j = ns; j <= ne ; j++) {
			if (j != i) {
				th *= (ah - x[j]) / (x[i] - x[j]);
			}
		}
		h += th;
	}
	return h / 60.0;
}

// 距離と標高差から山等のおよその仰角を求める(距離 dist 標高差 dh いずれも m
// 鈴木敬信著、地文航法より（気差と地球の曲率を考慮した略算式）
// 戻り値は仰角（radian)
function GetMountainElevation(dist,dh) {
	return Math.atan2((dh / dist - 11 * dist / 26 / 6.370e6),1.0);
}

// -->