matlab-convert-aerospace-coordinates

航空宇宙の基礎

単位変換、時刻変換、座標変換、回転表現を行う Aerospace Toolbox の主要な関数です。

使用する場面

• 航空宇宙分野の単位系（SI、ヤード・ポンド法、海里単位系）の間で変換する
• ユリウス日、修正ユリウス日、10進年、TDB 時刻を計算する
• ECEF、ECI、LLA、または平面地球座標系の間で位置を変換する
• 測地緯度と地心緯度の間で変換する
• 座標変換（ECI↔ECEF、ECEF↔NED、機体↔風、機体↔安定座標系）用の DCM を構築する
• オイラー角、DCM、クォータニオン、ロドリゲスベクトルの間で変換する
• 回転計算、補間、合成に quaternion オブジェクトを使用する
• データがどの座標系で表現されているかを特定する

使用しない場面

• 大気モデル（atmosisa、atmoscoesa）— matlab-compute-aerospace-environment を使用してください
• 対気速度の補正（correctairspeed）— ここでは扱いません
• 軌道伝播または衛星のマヌーバ — Aerospace Toolbox、Aerospace Blockset、または Satellite Communications Toolbox を使用してください
• ナビゲーション固有の変換（lla2enu、lla2ned）— Navigation Toolbox を使用してください
• Simulink ブロック — Aerospace Blockset とともに /model-based-design-core:building-simulink-models を使用してください

ワークフロー

1. 座標系を特定する — データがどの座標系にあり、どの座標系が必要かを判断します。座標系の定義と選択ガイドについては、references/coordinate-systems.md を参照してください。
2. 最初に単位を変換する — 変換関数を呼び出す前に、入力が関数で想定される単位と一致していることを確認します。
3. 変換を適用する — 適切な関数を正しい引数順序で使用します。
4. 検証する — 結果を元の座標系に戻して確認します。誤差は < 1e-10 になるはずです。

主要な関数

単位変換

関数	変換対象	単位
convlength	長さ	'ft'、'm'、'km'、'in'、'mi'、'naut mi'
convvel	速度	'ft/s'、'm/s'、'km/s'、'in/s'、'km/h'、'mph'、'kts'、'ft/min'
convang	角度	'deg'、'rad'、'rev'
convacc	加速度	'ft/s^2'、'm/s^2'、'km/s^2'、'in/s^2'、'km/h-s'、'mph/s'、'G''s'
convangacc	角加速度	'deg/s^2'、'rad/s^2'、'rpm/s'
convangvel	角速度	'deg/s'、'rad/s'、'rpm'
convforce	力	'N'、'lbf'
convmass	質量	'kg'、'lbm'、'slug'
convpres	圧力	'Pa'、'psi'、'psf'、'atm'
convtemp	温度	'K'、'R'、'F'、'C'
convdensity	密度	'kg/m^3'、'slug/ft^3'、'lbm/ft^3'、'lbm/in^3'

すべての変換関数で、同じシグネチャを使用します：output = convXXX(value, fromUnit, toUnit)

時刻変換

関数	目的	対応開始リリース
juliandate	暦日 → ユリウス日	R2006b
mjuliandate	暦日 → 修正ユリウス日（JD − 2400000.5）	R2006b
decyear	暦日 → 10進年	R2006b
leapyear	うるう年かどうかを判定	R2006b
tdbjuliandate	地球時 → TDB ユリウス日	R2015a

座標変換

関数	変換元	変換先	対応開始
lla2ecef	LLA（測地）	ECEF	R2006b
ecef2lla	ECEF	LLA（測地）	R2006b
lla2eci	LLA	ECI	R2014a
eci2lla	ECI	LLA	R2014a
ecef2eci	ECEF（位置／速度／加速度）	ECI	R2019a
eci2ecef	ECI（位置／速度／加速度）	ECEF	R2019a
eci2aer	ECI	AER（方位角、仰角、距離）	R2015a
lla2flat	LLA	平面地球	R2011a
flat2lla	平面地球	LLA	R2011a
geod2geoc	測地緯度	地心緯度	R2006b
geoc2geod	地心緯度	測地緯度	R2006b
ijk2keplerian	ECI 位置／速度	ケプラー要素	R2019a
keplerian2ijk	ケプラー要素	ECI 位置／速度	R2019a
siderealTime	UTC ユリウス日	GMST（度）／GAST（秒）	R2021a
dcmeci2ecef	—	ECI から ECEF への DCM	R2013b
dcmecef2ned	—	ECEF から NED への DCM	R2006b
dcm2latlon	ECEF から NED への DCM	緯度／経度	R2006b
dcmbody2wind	アルファ、ベータ	機体から風への DCM	R2006b
dcm2alphabeta	機体から風への DCM	アルファ、ベータ	R2006b
dcmbody2stability	アルファ	機体から安定軸への DCM	R2022a

回転表現

関数	変換元	変換先	対応開始
angle2dcm	オイラー角	DCM	R2006b
dcm2angle	DCM	オイラー角	R2006b
angle2quat	オイラー角	クォータニオン（1×4）	R2006b
quat2angle	クォータニオン（1×4）	オイラー角	R2007b
dcm2quat	DCM	クォータニオン（1×4）	R2006b
quat2dcm	クォータニオン（1×4）	DCM	R2006b
angle2rod	オイラー角	ロドリゲスベクトル	R2017a
rod2angle	ロドリゲスベクトル	オイラー角	R2017a
dcm2rod	DCM	ロドリゲスベクトル	R2017a
rod2dcm	ロドリゲスベクトル	DCM	R2017a
quat2rod	クォータニオン（1×4）	ロドリゲスベクトル	R2017a
rod2quat	ロドリゲスベクトル	クォータニオン（1×4）	R2017a

クォータニオンオブジェクト

メソッド	用途
quaternion(E,'eulerd',RS,PF)	オイラー角（度）から作成
quaternion(E,'euler',RS,PF)	オイラー角（ラジアン）から作成
quaternion(RM,'rotmat',PF)	回転行列から作成
quaternion(RV,'rotvec')	回転ベクトル（ラジアン）から作成
compact(q)	[w x y z] 配列を抽出
eulerd(q,RS,PF)	オイラー角（度）に変換
euler(q,RS,PF)	オイラー角（ラジアン）に変換
rotmat(q,PF)	回転行列に変換
rotvec(q) / rotvecd(q)	回転ベクトル（ラジアン／度）に変換
rotatepoint(q,pts)	点を回転（能動回転）
rotateframe(q,pts)	座標系を回転（受動回転）
normalize(q)	単位クォータニオンに正規化
slerp(q1,q2,t)	球面線形補間
meanrot(q)	配列の平均回転
dist(q1,q2)	角距離（ラジアン）
angvel(q,dt,PF)	クォータニオン配列から角速度を計算
randrot(n)	一様分布のランダム回転

クォータニオン演算（配列ベース）

関数	用途	対応開始
quatmultiply(q,r)	クォータニオン積（回転を合成）	R2006b
quatconj(q)	共役（ベクトル部の符号を反転）	R2006b
quatinv(q)	逆元（共役／ノルム²）	R2006b
quatnormalize(q)	単位クォータニオンに正規化	R2006b
quatnorm(q)	ノルムの二乗（q·q）	R2006b
quatmod(q)	絶対値（ノルムの平方根）	R2006b
quatrotate(q,v)	クォータニオンでベクトルを回転	R2006b
quatdivide(q,r)	クォータニオンをクォータニオンで除算	R2006b
quatinterp(p,q,f,method)	補間（'slerp'、'lerp'、'nlerp'）	R2016a
quatexp(q)	クォータニオンの指数関数	R2016a
quatlog(q)	クォータニオンの自然対数	R2016a
quatpower(q,pow)	クォータニオンをべき乗	R2016a

パターン

単位変換

% Always: convXXX(value, 'from', 'to')
alt_m = convlength(35000, 'ft', 'm');        % 10668.0 m
speed_ms = convvel(250, 'kts', 'm/s');       % 128.61 m/s
angle_rad = convang(45, 'deg', 'rad');       % 0.7854 rad
accel_g = convacc(9.81, 'm/s^2', 'G''s');   % 1.0 G
omega_rpm = convangvel(360, 'deg/s', 'rpm'); % 60 rpm

時刻変換

% Julian Date from components (year, month, day, hour, min, sec)
jd = juliandate(2024, 6, 15, 12, 0, 0);     % 2460477.0

% Modified Julian Date
mjd = mjuliandate(2024, 6, 15, 12, 0, 0);   % 60476.5

% From datetime objects
dt = datetime(2024, 6, 15, 12, 0, 0);
jd = juliandate(dt);

% Decimal year
dy = decyear(2024, 6, 15);                   % 2024.4536

% TDB Julian Date from Terrestrial Time [yr mo day hr min sec]
tt = [2024 6 15 12 0 0];
jdTDB = tdbjuliandate(tt);

LLA ↔ ECEF

% LLA is [latitude_deg, longitude_deg, altitude_m]
lla = [40, -74, 0];
ecef = lla2ecef(lla);           % [1348613.0, -4703172.4, 4077985.6] m

% Round-trip verification
lla_check = ecef2lla(ecef);     % [40.0, -74.0, 0.0]

LLA ↔ ECI（時間依存）

% ECI transforms require UTC time
lla = [40, -74, 1000];
utc = [2024 6 15 12 0 0];
posECI = lla2eci(lla, utc);

% Back to LLA
lla_check = eci2lla(posECI, utc);

ケプラー軌道要素 ↔ ECI（IJK）

% Convert ECI position/velocity to Keplerian orbital elements
rijk = [-2981784; 5207055; 3161595];   % ECI position (m), 3x1
vijk = [-3384; -4887; 4843];           % ECI velocity (m/s), 3x1
[a, ecc, incl, RAAN, argp, nu] = ijk2keplerian(rijk, vijk);
% a=m, angles in degrees

% Convert Keplerian elements back to ECI position/velocity
[rijk_back, vijk_back] = keplerian2ijk(a, ecc, incl, RAAN, argp, nu);

% For a different central body (e.g., Mars)
[a_mars, ecc_mars, incl_mars, RAAN_mars, argp_mars, nu_mars] = ...
    ijk2keplerian(rijk, vijk, CentralBody="Mars");

恒星時

% Compute Greenwich Mean Sidereal Time (GMST) and Apparent (GAST)
utcJD = juliandate(2019, 1, 4, 12, 0, 0);  % UTC as Julian date
dUT1 = 0;    % UT1-UTC offset (seconds)
dAT  = 37;   % TAI-UTC offset (leap seconds)
[thGMST, thGAST] = siderealTime(utcJD, dUT1, dAT);
% thGMST in degrees, thGAST in seconds

平面地球近似

% Good for short-range simulations (< ~100 km from reference)
llo = [40, -74];    % reference lat/lon (deg)
psio = 0;           % angular direction of flat Earth x-axis (rad, 0=North)
href = 0;           % reference height (m)

lla_point = [40.01, -73.99, 100];
flatPos = lla2flat(lla_point, llo, psio, href);  % [x, y, z] in meters
lla_back = flat2lla(flatPos, llo, psio, href);

オイラー角 ↔ DCM

% CRITICAL: angle2dcm expects RADIANS, not degrees
yaw = deg2rad(30); pitch = deg2rad(10); roll = deg2rad(5);
dcm = angle2dcm(yaw, pitch, roll, 'ZYX');

% Extract angles back (returns radians)
[y, p, r] = dcm2angle(dcm, 'ZYX');
fprintf('Yaw=%.1f, Pitch=%.1f, Roll=%.1f deg\n', rad2deg(y), rad2deg(p), rad2deg(r));

DCM ↔ クォータニオン（Aerospace Toolbox 形式）

% Aerospace Toolbox quaternion format: [q0 q1 q2 q3] = [scalar, vector]
q = dcm2quat(dcm);       % 1x4, scalar-first
dcm_back = quat2dcm(q);  % 3x3

% Verify orthogonality
err = max(abs(dcm' * dcm - eye(3)), [], 'all');
assert(err < 1e-14, 'DCM is not orthogonal');

クォータニオンオブジェクト（最新のアプローチ）

% Create from Euler angles — specify sequence AND frame type
q = quaternion([30 10 5], 'eulerd', 'ZYX', 'frame');

% Rotate a point (active rotation)
pt = [1 0 0];
pt_rotated = rotatepoint(q, pt);   % [0.8529, 0.4924, -0.1736]

% Compose rotations by multiplication
q_total = q2 * q1;  % applies q1 first, then q2

% Always normalize after arithmetic accumulation
q = normalize(q);

クォータニオン補間（SLERP）

% Smooth interpolation between two orientations
q1 = quaternion([0 0 0], 'eulerd', 'ZYX', 'frame');
q2 = quaternion([90 0 0], 'eulerd', 'ZYX', 'frame');
t = linspace(0, 1, 5)';
q_interp = slerp(q1, q2, t);

% Verify: yaw progresses linearly for pure yaw rotation
e = eulerd(q_interp, 'ZYX', 'frame');
% e(:,1) = [0, 22.5, 45, 67.5, 90]

クォータニオン演算（配列ベースの関数）

% Compose two rotations: q2 applied after q1
q1 = angle2quat(deg2rad(30), 0, 0, 'ZYX');
q2 = angle2quat(0, deg2rad(10), 0, 'ZYX');
q_total = quatmultiply(q2, q1);  % q2 * q1 (apply q1 first)

% Inverse rotation
q_inv = quatinv(q1);
q_identity = quatmultiply(q1, q_inv);  % [1 0 0 0]

% Rotate a vector
v = [1 0 0];
v_rot = quatrotate(q1, v);  % rotate v by q1

% Normalize after accumulation
q_accumulated = quatnormalize(q_total);

% Interpolation (slerp, lerp, nlerp)
q_mid = quatinterp(q1, q2, 0.5, 'slerp');

% Relative rotation: what rotation takes q1 to q2?
q_rel = quatdivide(q2, q1);  % q_rel such that q2 = q_rel * q1

% Exponential/logarithm (useful for angular velocity integration)
q_log = quatlog(q1);         % maps to tangent space
q_back = quatexp(q_log);     % back to quaternion

% Fractional rotation (half the rotation of q1)
q_half = quatpower(q1, 0.5);

ロドリゲスベクトル

% Compact 3-element representation (singular at 180 deg)
rod = dcm2rod(dcm);           % 1x3 vector
dcm_back = rod2dcm(rod);

% Convert between all representations
rod = quat2rod(q_array);      % q is M×4 (scalar-first)
q_back = rod2quat(rod);

機体座標系から風軸座標系および安定軸座標系への変換

% Body-to-wind DCM from angle of attack and sideslip
alpha = deg2rad(5);
beta = deg2rad(2);
dcm_bw = dcmbody2wind(alpha, beta);

% Extract alpha/beta from a DCM
[alpha_out, beta_out] = dcm2alphabeta(dcm_bw);

% Body-to-stability (alpha only, no sideslip)
dcm_bs = dcmbody2stability(alpha);

ECEF から NED への座標系 DCM

% Get the DCM to rotate vectors from ECEF to local NED
lat = deg2rad(40);
lon = deg2rad(-74);
dcm_ecef2ned = dcmecef2ned(lat, lon);

% Transform an ECEF velocity to NED
v_ecef = [10; 20; 30];
v_ned = dcm_ecef2ned * v_ecef;

ECI から AER（方位角、仰角、斜距離）への変換

% Compute look angles from a ground station to a satellite in ECI
posECI = [-2981784, 5207055, 3161595];  % satellite ECI position (m)
utc = [2019 1 4 12 0 0];               % observation time
lla0 = [28.5, -80.5, 0];               % ground station [lat, lon, alt] (deg, deg, m)

aer = eci2aer(posECI, utc, lla0);
% aer = [azimuth_deg, elevation_deg, slant_range_m]

規約

• angle2dcm と angle2quat はラジアンを想定します — 事前に必ず deg2rad() または convang で変換してください
• 回転シーケンスの既定値は 'ZYX'（ヨー・ピッチ・ロール）です — 明確にするため、常に明示的に指定してください
• Aerospace Toolbox のクォータニオン形式は [scalar, i, j, k] です — 1×4 配列関数はスカラーを先頭にします
• quaternion オブジェクトもスカラーを先頭にします（compact は [w x y z] を返します）
• rotatepoint と rotateframe の違い — rotatepoint は点を回転させ（能動回転）、rotateframe は座標系を回転させます（受動回転）。この2つは互いに逆の変換です。
• LLA の順序は [lat, lon, alt] です。単位は度とメートルです
• geod2geoc と geoc2geod は緯度にラジアンを、標高／半径にメートルを想定します
• ECI 関数には UTC 時刻が必要です — 結果は地球の自転に応じて変化します
• 演算後はクォータニオンを正規化してください。ドリフトを防止できます
• ロドリゲスベクトルは 180° で特異点があります — 任意の回転にはクォータニオンを使用してください

よくある間違い

間違い	問題となる理由	正しい方法
angle2dcm に度を渡す	関数はラジアンを想定します	deg2rad() または convang(val,'deg','rad') を使用します
回転シーケンスを省略する	既定の ZYX が使用する規約と一致しない可能性があります	常に 'ZYX'、'XYZ' などを明示的に渡します
クォータニオンが [x y z w] だと思い込む	Aerospace Toolbox は [w x y z]（スカラーが先頭）を使用します	形式を確認し、compact(q) で確かめます
rotateframe のつもりで rotatepoint を使う	互いに逆の変換なので、誤った方を使うと回転が反転します	能動回転には rotatepoint、受動回転には rotateframe を使用します
クォータニオン演算後に正規化しない	クォータニオンのドリフトでノルムが 1 でなくなり、回転が歪みます	回転を累積した後に normalize(q) を呼び出します
geod2geoc に度を渡す	緯度にはラジアンが必要です	geod2geoc(deg2rad(lat), alt) のように変換します
ECI 変換で UTC を指定し忘れる	ECI 位置はその時点での地球の自転に依存します	常に [yr mo day hr min sec] を指定します
180° 付近の回転にロドリゲスベクトルを使う	ロドリゲスベクトルは π で特異点があります	代わりにクォータニオン表現を使用します
ピッチが ±90° 付近のジンバルロックを無視する	オイラー角の抽出で自由度が 1 つ失われます	計算にはクォータニオンまたは DCM を直接使用します

ジンバルロック

2番目の回転が ±90° に達すると、オイラー角表現では自由度が 1 つ失われます（ZYX の場合、ピッチ = ±90°）。次のような症状が現れます。

• ピッチが ±90° に近いとき、dcm2angle が予期しないヨー／ロールの値を返す
• 特異点付近でオイラー角を補間すると、経路が不規則になる

解決策： 計算と補間にはクォータニオンを使用します。オイラー角への変換は、表示または人間が解釈する場合に限って行ってください。

% Detect gimbal lock risk
[~, pitch, ~] = dcm2angle(dcm, 'ZYX');
if abs(pitch) > deg2rad(85)
    warning('Near gimbal lock — use quaternion representation');
end

参考資料

• 座標系の詳細な定義（ECEF、ECI、NED、機体、風軸、安定軸）、軸の規約、および適切な座標系を選択するためのガイドについては、references/coordinate-systems.md を参照してください。
• 回転に関する詳細な数学（DCM の特性、クォータニオン代数、ジンバルロックの理論、および各種表現間の変換経路）については、references/rotation-conventions.md を参照してください。
• 単位文字列の完全な一覧と、連鎖変換については、references/unit-conversion-guide.md を参照してください。

---

Copyright 2026 The MathWorks, Inc.

---