matlab-compute-aerospace-environment

航空宇宙環境

Aerospace Toolbox の関数を使用して、航空宇宙機の解析に必要な環境特性（大気、重力、風、磁場、ジオイド、宇宙天気、惑星暦、地球姿勢）を計算します。

使用する場面

• 大気特性（温度、圧力、密度、音速）の計算
• 任意の惑星の重力ベクトルまたは加速度の計算
• 高度における水平風のモデル化
• 航法やコンパス補正に必要な磁場成分の取得
• ジオイド高または地心半径の計算
• NRLMSISE-00 の入力に使う宇宙天気データの読み込み
• 惑星または月の位置（惑星暦）の計算
• 地球姿勢パラメーター（極運動、章動、UT1-UTC、CIP）の取得
• 大気モデル、ISA、COESA、NRLMSISE、圧力、温度、密度、音速、重力モデル、WGS84、風モデル、HWM、磁場モデル、WMM、IGRF、ジオイド、宇宙天気、太陽フラックス、F10.7、Ap 指数、惑星暦、月の位置、地球の章動、極運動、UT1、IERS に言及するプロンプト

使用しない場面

• 座標系の変換または回転 — matlab-convert-aerospace-coordinates を使用してください
• 軌道力学または軌道伝播（位置には惑星暦を使用し、軌道伝播には使用しないでください）
• 空力係数の計算
• Simulink の環境モデルブロック — Aerospace Blockset を使用してください

ワークフロー

1. 必要な量を特定する — 以下の判断表を使用してください
2. 関数を呼び出す — このスキルのパターンに従って、正しい構文で呼び出してください
3. 結果を確認する — 単位と値の大きさが物理的に妥当であることを確認してください

判断表

必要な量	関数	主な入力
標準大気（簡易）	atmosisa	高度（m）
1976 COESA 大気	atmoscoesa	高度（m）
NRLMSISE-00（詳細な高層大気）	atmosnrlmsise00	高度、緯度、経度、年、通日、UT秒
非標準大気（MIL-STD）	atmosnonstd	高度 + 位置引数の文字列
CIRA 1986 基準大気	atmoscira	緯度、ctype、coord、月
温度減率大気（カスタム）	atmoslapse	高度（m）+ 9 個の物理パラメーター
気圧高度	atmospalt	圧力（Pa）
水平風（HWM07/14）	atmoshwm	緯度、経度、高度 + 名前と値
球面調和重力（任意の惑星）	gravitysphericalharmonic	PCPF [x,y,z]（m）
WGS84 重力（地球、測地）	gravitywgs84	h、緯度（+ 経度、method、flags）
帯状調和重力（任意の惑星）	gravityzonal	PCPF [x,y,z]（m）
遠心加速度	gravitycentrifugal	PCPF [x,y,z]（m）
WMM 磁場	wrldmagm	高度（m）、緯度、経度、10 進年
IGRF 磁場	igrfmagm	高度（m）、緯度、経度、10 進年、世代
ジオイド高（起伏）	geoidheight	緯度、経度
地心半径	geocradius	地心緯度（度）
宇宙天気 CSV の読み込み	aeroReadSpaceWeatherData	CSV ファイルのパス
太陽フラックス / Ap 指数の抽出	fluxSolarAndGeomagnetic	datetime または [年、通日、UTC秒]、MAT ファイル
惑星 / 月の位置と速度	planetEphemeris	ユリウス日、中心、対象
地球の章動角	earthNutation	ユリウス日
月の秤動角	moonLibration	ユリウス日
地球の極運動	polarMotion	UTC（ユリウス日）
天球中間極の補正	deltaCIP	UTC（ユリウス日）
UT1 と UTC の差	deltaUT1	UTC（ユリウス日）
IERS 地球姿勢データの読み込み	aeroReadIERSData	フォルダーのパス

使用パターン

標準大気および COESA 大気

% International Standard Atmosphere
[T, a, P, rho] = atmosisa(1000);

% 1976 COESA (valid 0-1000 km)
[T, a, P, rho] = atmoscoesa(1000);

気圧高度（atmospalt）

国際標準大気を使用して、圧力（Pa）を高度（m）に変換します。

% Pressure altitude at standard sea-level pressure
alt = atmospalt(101325);  % returns 0 m

% Pressure altitude at multiple pressures
alt = atmospalt([101325, 79501, 54048, 26500]);

% Typical use: convert measured pressure to altitude
measuredPressure_Pa = 75000;
pressureAltitude_m = atmospalt(measuredPressure_Pa);

入力：圧力（パスカル）。出力：1976 COESA に基づく幾何高度（メートル）。

非標準大気 (atmosnonstd)

位置文字列引数を使用します — 名前と値のペアではありません。

プロファイルタイプ (単一の高度極値):

[T, a, P, rho] = atmosnonstd(height, 'Profile', extremeParam, frequency, extremeAltitude)

エンベロープタイプ (高度範囲の極値 — extremeAltitude 引数はありません):

[T, a, P, rho] = atmosnonstd(height, 'Envelope', extremeParam, frequency)

省略可能な末尾の引数: action ('Warning'|'Error'|'None')、specification ('310'|'210c')。

% Profile: high density, 1% of time, at 5 km altitude (extremeAltitude is NUMERIC)
[T, a, P, rho] = atmosnonstd(5000, 'Profile', 'High density', '1%', 5);

% Envelope: high pressure, 20% of time, MIL-STD-210C
[T, a, P, rho] = atmosnonstd([1000; 11000; 20000], 'Envelope', ...
    'High pressure', '20%', 'None', '210c');

有効な extremeParam: 'High temperature'、'Low temperature'、'High density'、'Low density'、'High pressure'、'Low pressure'

有効な frequency: 'Extreme values'、'1%'、'5%'、'10%'、'20%'

有効な extremeAltitude (プロファイルのみ、数値): 5、10、20、30、40

CIRA 1986 モデル

% By geopotential height, monthly mean, October
[T, P, zonalWind] = atmoscira(45, 'GPHeight', 20000, 'Monthly', 10);

% By pressure level
[T, alt, zonalWind] = atmoscira(45, 'Pressure', 101300, 'Monthly', 1);

引数: (latitude, ctype, coord, mtype, month)。ここで ctype は 'Pressure' または 'GPHeight' です。

NRLMSISE-00

% Basic call (uses default flux values)
[T, rho] = atmosnrlmsise00(altitude, latitude, longitude, year, dayOfYear, UTseconds);

% With flux data and no anomalous oxygen
[T, rho] = atmosnrlmsise00(altitude, lat, lon, year, dayOfYear, UTsec, ...
    f107Average, f107Daily, magneticIndex, 'NoOxygen');

出力: T は N×2 [外気圏温度、局所温度] です。rho は N×9 (密度) です。

気温減率大気 (atmoslapse)

10 個の位置引数を必要とします — 高度だけではありません。すべての物理パラメータを明示的に指定する必要があります。

% TEMPLATE — not executable (shows calling convention)
[T, a, P, rho] = atmoslapse(altitude, g, gamma, R, lapseRate, ...
    heightTroposphere, heightTropopause, density0, pressure0, temperature0);

% Example: ISA conditions
[T, a, P, rho] = atmoslapse([0 5000 11000 20000], ...
    9.80665, 1.4, 287.0531, 0.0065, 11000, 20000, 1.225, 101325, 288.15);

引数: (height, g, heatRatio, gasConstant, lapseRate, hTroposphere, hTropopause, rho0, P0, T0)。カスタムの惑星大気や非標準の海面条件に使用します。

水平風モデル (atmoshwm)

名前と値の引数を使用します。M×2 配列 [子午線方向、帯状方向] を m/s で返します。

% Total wind with Ap index
wind = atmoshwm(-45, -85, 25000, model='total', day=150, ...
    seconds=11*3600, apindex=80);

% Quiet wind at multiple altitudes (day/seconds must match array size)
wind = atmoshwm([50; 50], [-20; -20], [100000; 150000], ...
    model='quiet', day=[30; 30], seconds=[0; 0]);

% Disturbed wind
wind = atmoshwm(70, -65, 150000, model='disturbance', day=166, seconds=0);

名前と値のオプション: model ('quiet'|'disturbance'|'total')、day、seconds、apindex、version ('14'|'07')。

重力: 球面調和関数

入力は惑星中心惑星固定 (PCPF) [x,y,z] 座標 (メートル単位) です — 緯度/経度/高度ではありません。

% Earth: 25,000 m over South Pole (z-axis = polar axis)
[gx, gy, gz] = gravitysphericalharmonic([0, 0, -(6356752 + 25000)]);

% Mars GMM2B: equator + pole
pos = [3396200 + 15000, 0, 0;   % equator
       0, 0, 3376200 + 11000];  % north pole
[gx, gy, gz] = gravitysphericalharmonic(pos, 'GMM2B');

モデル: 'EGM2008' (デフォルト)、'EGM96'、'LP100K'、'LP165P'、'GMM2B'、'EIGENGL04C'、'Custom'。

注: 扁球体上の極位置では「Radial position is less than equatorial radius」という警告が発生します — これは想定される動作であり (極半径 < 赤道半径)、結果は依然として有効です。

重力: WGS84

測地座標を使用します (高度はメートル、緯度/経度は度)。

オプションの効果を伴う Exact 法では、フラグを単一のベクトル [noatm, nocent, prec, jd] として渡します:

% Taylor series (default) — only h and lat, no method string
g = gravitywgs84(1000, 45);

% CloseApprox — requires lon as 3rd argument
g = gravitywgs84(1000, 45, 20, 'CloseApprox');

% Exact with precession — requires lon + flag vector
[gn, gt] = gravitywgs84(1000, 0, 20, 'Exact', [false, false, true, 2451545]);

呼び出し規約:

• デフォルト (テイラー級数): gravitywgs84(h, lat) — 'TaylorSeries' を文字列として渡さないでください
• CloseApprox: gravitywgs84(h, lat, lon, 'CloseApprox') — lon が必要です
• Exact: gravitywgs84(h, lat, lon, 'Exact', [noatm, nocent, prec, jd]) — lon + フラグベクトルが必要です

フラグベクトル [noatm, nocent, prec, jd]:

• noatm: true = 大気補正を除外
• nocent: true = 遠心効果を除外
• prec: true = 歳差を含める
• jd: ユリウス日 (0 = 歳差の日付なし)

重力: 帯球関数

入力はPCPF [x,y,z] (メートル単位) です。次数は数値スカラーです。

% Earth, 4th degree, at equator surface
Re = 6378137;
[gx, gy, gz] = gravityzonal([Re, 0, 0], 'Earth', 4);

% Mars at two positions
[gx, gy, gz] = gravityzonal([3396200+15000, 0, 0; 0, 0, 3376200+11000], 'Mars');

重力: 遠心力

入力はPCPF [x,y,z] (メートル単位) です。

% Earth at equator surface
[gx, gy, gz] = gravitycentrifugal([6378137, 0, 0]);

% Mars
[gx, gy, gz] = gravitycentrifugal([3396200, 0, 0], 'Mars');

磁場モデル

wrldmagm と igrfmagm はどちらも高度をメートル、緯度/経度を度、および小数年で受け取ります。暦日を変換するには decyear を使用します。

% WMM-2025
[XYZ, H, D, I, F] = wrldmagm(1000, 40.44, -79.99, decyear(2025,7,4), '2025');

% IGRF-14
[XYZ, H, D, I, F] = igrfmagm(3000, 41.50, -81.69, decyear(2028,7,4), 14);

出力: XYZ = [北、東、下] (nT); H = 水平強度; D = 偏角 (度); I = 傾角 (度); F = 全強度 (nT)。

ジオイド高と地心半径

% Geoid undulation (EGM96 default) — longitude must be [0, 360]
N = geoidheight(42.4, 289.0);

% EGM2008
N = geoidheight(42.4, 289.0, 'EGM2008');

% Geocentric radius at multiple latitudes (input: geocentric lat in degrees)
r = geocradius([0, 45, 90]);           % default WGS84
r = geocradius([0, 45, 90], 'WGS84');  % explicit model

惑星暦と地球姿勢

暦関数は aeroDataPackage からのデータを必要とします。すべて時間入力としてユリウス日を受け取ります。

% Planet/Moon position and velocity (km, km/s by default)
jd = juliandate(2022, 10, 27, 12, 0, 0);
[pos, vel] = planetEphemeris(jd, 'Earth', 'Moon');
[pos, vel] = planetEphemeris(jd, 'SolarSystem', 'Mars', '430', 'AU');

% Earth nutation (radians)
[angles, rates] = earthNutation(jd);  % angles = [dpsi, depsilon]

% Moon libration (radians)
[angles, rates] = moonLibration(jd);  % 3 Euler angles + rates

% Earth polar motion (radians) — input is UTC as Julian date
utc = juliandate(2020, 1, 1);
[pm, pmError] = polarMotion(utc);  % pm = [xp, yp]

% Celestial Intermediate Pole adjustment (radians)
[dcip, dcipError] = deltaCIP(utc);  % dcip = [dX, dY]

% UT1-UTC difference (seconds)
[dut1, dut1Error] = deltaUT1(utc);

有効な planetEphemeris のターゲット/中心: 'Sun'、'Mercury'、'Venus'、'Earth'、'Moon'、'Mars'、'Jupiter'、'Saturn'、'Uranus'、'Neptune'、'Pluto'、'SolarSystem'、'EarthMoon'。

モデル: '405' (デフォルト)、'421'、'423'、'430'、'432t'。

注: polarMotion、deltaCIP、deltaUT1 は aeroiersdata.mat (同梱) またはカスタムソースからの IERS データを使用します。更新された IERS データを取得するには:

% Download latest Earth orientation data and save as MAT file
iersFile = aeroReadIERSData(pwd);
% Then pass to Earth orientation functions via Source name-value
[pm, pmError] = polarMotion(utc, Source=iersFile);
[dcip, dcipError] = deltaCIP(utc, Source=iersFile);
[dut1, dut1Error] = deltaUT1(utc, Source=iersFile);

宇宙天気パイプライン

現実的な NRLMSISE-00 入力値を使用する完全なパイプライン：

% Step 1: Read CSV (downloaded from CelesTrack) into MAT file
% TEMPLATE — requires user-provided CSV file
matFile = aeroReadSpaceWeatherData('spaceWeatherFile.csv');

% Step 2: Extract flux and magnetic indices for a specific date/time
dt = datetime(2022, 10, 27, 6, 30, 0, 'TimeZone', 'UTC');
[f107avg, f107daily, magIdx] = fluxSolarAndGeomagnetic(dt, matFile, ...
    F107ExtrapMethod='constant', MagneticIndexExtrapMethod='constant');

% Step 3: Feed into NRLMSISE-00
dayOfYear = day(dt, 'dayofyear');
UTsec = seconds(timeofday(dt));
[T, rho] = atmosnrlmsise00(400e3, 45, -50, dt.Year, dayOfYear, UTsec, ...
    f107avg, f107daily, magIdx, 'NoOxygen');

注意：

• aeroReadSpaceWeatherData はユーザーが用意した CSV（CelesTrack 形式）を受け取ります。ツールボックスにはサンプルは付属していません
• fluxSolarAndGeomagnetic で使用する宇宙天気データを含む MAT ファイルを返します
• CSV が指定されていない場合は、事前構築済みの aeroSpaceWeatherData.mat（which('aeroSpaceWeatherData.mat') 経由）を利用できます
• fluxSolarAndGeomagnetic は (datetime, matFile) または (year, dayOfYear, UTCseconds, matFile) を受け付けます

よくある間違い

間違い	失敗する理由	正しい方法
gravitysphericalharmonic に緯度／経度／高度を渡す	PCPF の [x,y,z]（単位：メートル）が必要です	直交座標に変換します：赤道面 = [Re, 0, 0]、極 = [0, 0, Rp]
atmosnonstd に名前と値のペアを渡す	名前と値ではなく、位置指定文字列を使用します	atmosnonstd(h, 'Envelope', 'High pressure', '20%')
gravitywgs84 Exact に個別のブール値引数を渡す	フラグは単一のベクトルで指定する必要があります	gravitywgs84(h, lat, lon, 'Exact', [noatm, nocent, prec, jd])
メソッド文字列として 'TaylorSeries' を渡す	メソッド引数が必要な場合、3 番目の引数です	2 引数の形式を使用します：gravitywgs84(h, lat) — メソッド文字列は不要です
経度を指定せずに gravitywgs84(h, lat, 'CloseApprox') を呼び出す	CloseApprox には経度が必要です	gravitywgs84(h, lat, lon, 'CloseApprox')
Envelope 型に extremeAltitude を含める	Profile 型でのみ有効です	Envelope では省略します：atmosnonstd(h, 'Envelope', param, freq)
風に atmosnrlmsise00 を使用する	関数が違います	水平風には atmoshwm を使用します
遠心力の影響を得るために gravitywgs84 を呼び出す	遠心力のみを計算するモードはありません	gravitycentrifugal([x, y, z]) を使用します
関数名 gravityzonalharmonic を誤って使用する	存在しません	gravityzonal を使用します
度数を文字列で gravityzonal に渡す	数値のスカラーが必要です	'4thDegree' ではなく gravityzonal(coords, 'Earth', 4) を使用します
aeroReadSpaceWeatherData の出力を無視する	宇宙天気データを含む MAT ファイルを返します	返された MAT ファイルを fluxSolarAndGeomagnetic に渡します

規則

• 高度の単位： 大気、重力、磁気のすべての関数でメートルを使用します（km やフィートではありません）
• 緯度／経度： 測地関数では度数を使用します。gravitysphericalharmonic、gravityzonal、gravitycentrifugal では PCPF の [x,y,z]（単位：メートル）を使用します
• 10 進年： 磁気モデルの日付入力には decyear(year, month, day) を使用します
• 時刻： day = 年通算日（1～366）；seconds = UTC の午前 0 時からの秒数
• PCPF 座標： x 軸は 0N/0E を通り、y 軸は 0N/90E を通り、z 軸は北極を通ります
• geoidheight の経度： [0, 360] 度（東が正）で指定する必要があります。[-180, 180] から変換するには mod(lon, 360) を使用します
• PCPF の南極： [0, 0, -(Rp + altitude)]。Rp は極半径です
• PCPF の赤道： [Re + altitude, 0, 0]。Re は赤道半径です

---

Copyright 2026 The MathWorks, Inc.

---