From 7f847e027453810049cb2d87a62b91d8a335a9cd Mon Sep 17 00:00:00 2001 From: Claude Date: Wed, 19 Nov 2025 03:03:56 +0000 Subject: [PATCH] feat: add advanced STK features to OOTK Add comprehensive satellite analysis features inspired by STK: - Coverage Analysis: ground track calculation, area coverage grids, access statistics - Conjunction/Collision Analysis: close approach detection, probability of collision calculation using Chan, Foster, and Patera methods - Link Budget Analysis: RF link calculations including path loss, SNR, CNR, Eb/N0, atmospheric losses - Walker Constellations: Walker Delta and Star patterns with preset configurations (GPS, Galileo, Iridium, OneWeb) - Attitude Dynamics: quaternion-based attitude representation, Euler angles conversion, attitude propagation with torque models These additions significantly expand OOTK's capabilities for mission planning, constellation design, and spacecraft operations analysis. --- package-lock.json | 230 ++-------------- package.json | 1 + src/attitude/Attitude.ts | 281 +++++++++++++++++++ src/attitude/AttitudePropagator.ts | 233 ++++++++++++++++ src/attitude/index.ts | 21 ++ src/conjunction/ConjunctionAnalysis.ts | 255 ++++++++++++++++++ src/conjunction/ProbabilityOfCollision.ts | 255 ++++++++++++++++++ src/conjunction/index.ts | 21 ++ src/constellation/WalkerConstellation.ts | 304 +++++++++++++++++++++ src/constellation/index.ts | 20 ++ src/coverage/AccessStatistics.ts | 231 ++++++++++++++++ src/coverage/AreaCoverage.ts | 279 +++++++++++++++++++ src/coverage/GroundTrack.ts | 210 +++++++++++++++ src/coverage/index.ts | 22 ++ src/link_budget/LinkBudget.ts | 314 ++++++++++++++++++++++ src/link_budget/index.ts | 20 ++ src/main.ts | 10 + 17 files changed, 2495 insertions(+), 212 deletions(-) create mode 100644 src/attitude/Attitude.ts create mode 100644 src/attitude/AttitudePropagator.ts create mode 100644 src/attitude/index.ts create mode 100644 src/conjunction/ConjunctionAnalysis.ts create mode 100644 src/conjunction/ProbabilityOfCollision.ts create mode 100644 src/conjunction/index.ts create mode 100644 src/constellation/WalkerConstellation.ts create mode 100644 src/constellation/index.ts create mode 100644 src/coverage/AccessStatistics.ts create mode 100644 src/coverage/AreaCoverage.ts create mode 100644 src/coverage/GroundTrack.ts create mode 100644 src/coverage/index.ts create mode 100644 src/link_budget/LinkBudget.ts create mode 100644 src/link_budget/index.ts diff --git a/package-lock.json b/package-lock.json index a36a7437..4aee867d 100644 --- a/package-lock.json +++ b/package-lock.json @@ -14,6 +14,7 @@ "@babel/preset-env": "^7.26.9", "@babel/preset-typescript": "^7.26.0", "@types/jest": "^29.5.14", + "@types/node": "^24.10.1", "@typescript-eslint/eslint-plugin": "^7.18.0", "@typescript-eslint/parser": "^7.18.0", "auto-changelog": "^2.5.0", @@ -21,7 +22,6 @@ "babel-jest": "^28.1.3", "babel-loader": "^9.2.1", "eslint": "^8.57.1", - "eslint-plugin-jsdoc": "^48.11.0", "fix-esm-import-path": "^1.10.1", "jest": "^28.1.3", "ts-jest": "^28.0.8", @@ -1874,21 +1874,6 @@ "@jridgewell/sourcemap-codec": "^1.4.10" } }, - "node_modules/@es-joy/jsdoccomment": { - "version": "0.46.0", - "resolved": "https://registry.npmjs.org/@es-joy/jsdoccomment/-/jsdoccomment-0.46.0.tgz", - "integrity": "sha512-C3Axuq1xd/9VqFZpW4YAzOx5O9q/LP46uIQy/iNDpHG3fmPa6TBtvfglMCs3RBiBxAIi0Go97r8+jvTt55XMyQ==", - "dev": true, - "license": "MIT", - "dependencies": { - "comment-parser": "1.4.1", - "esquery": "^1.6.0", - "jsdoc-type-pratt-parser": "~4.0.0" - }, - "engines": { - "node": ">=16" - } - }, "node_modules/@esbuild/aix-ppc64": { "version": "0.25.1", "resolved": "https://registry.npmjs.org/@esbuild/aix-ppc64/-/aix-ppc64-0.25.1.tgz", @@ -3467,19 +3452,6 @@ "node": ">= 8" } }, - "node_modules/@pkgr/core": { - "version": "0.1.2", - "resolved": "https://registry.npmjs.org/@pkgr/core/-/core-0.1.2.tgz", - "integrity": "sha512-fdDH1LSGfZdTH2sxdpVMw31BanV28K/Gry0cVFxaNP77neJSkd82mM8ErPNYs9e+0O7SdHBLTDzDgwUuy18RnQ==", - "dev": true, - "license": "MIT", - "engines": { - "node": "^12.20.0 || ^14.18.0 || >=16.0.0" - }, - "funding": { - "url": "https://opencollective.com/unts" - } - }, "node_modules/@sinclair/typebox": { "version": "0.27.8", "resolved": "https://registry.npmjs.org/@sinclair/typebox/-/typebox-0.27.8.tgz", @@ -3657,10 +3629,14 @@ "dev": true }, "node_modules/@types/node": { - "version": "14.14.20", - "resolved": "https://registry.npmjs.org/@types/node/-/node-14.14.20.tgz", - "integrity": "sha512-Y93R97Ouif9JEOWPIUyU+eyIdyRqQR0I8Ez1dzku4hDx34NWh4HbtIc3WNzwB1Y9ULvNGeu5B8h8bVL5cAk4/A==", - "dev": true + "version": "24.10.1", + "resolved": "https://registry.npmjs.org/@types/node/-/node-24.10.1.tgz", + "integrity": "sha512-GNWcUTRBgIRJD5zj+Tq0fKOJ5XZajIiBroOF0yvj2bSU1WvNdYS/dn9UxwsujGW4JX06dnHyjV2y9rRaybH0iQ==", + "dev": true, + "license": "MIT", + "dependencies": { + "undici-types": "~7.16.0" + } }, "node_modules/@types/prettier": { "version": "2.7.3", @@ -4274,15 +4250,6 @@ "node": ">= 8" } }, - "node_modules/are-docs-informative": { - "version": "0.0.2", - "resolved": "https://registry.npmjs.org/are-docs-informative/-/are-docs-informative-0.0.2.tgz", - "integrity": "sha512-ixiS0nLNNG5jNQzgZJNoUpBKdo9yTYZMGJ+QgT2jmjR7G7+QHRCc4v6LQ3NgE7EBJq+o0ams3waJwkrlBom8Ig==", - "dev": true, - "engines": { - "node": ">=14" - } - }, "node_modules/arg": { "version": "4.1.3", "resolved": "https://registry.npmjs.org/arg/-/arg-4.1.3.tgz", @@ -4791,16 +4758,6 @@ "license": "MIT", "peer": true }, - "node_modules/comment-parser": { - "version": "1.4.1", - "resolved": "https://registry.npmjs.org/comment-parser/-/comment-parser-1.4.1.tgz", - "integrity": "sha512-buhp5kePrmda3vhc5B9t7pUQXAb2Tnd0qgpkIhPhkHXxJpiPJ11H0ZEU0oBpJ2QztSbzG/ZxMj/CHsYJqRHmyg==", - "dev": true, - "license": "MIT", - "engines": { - "node": ">= 12.0.0" - } - }, "node_modules/common-path-prefix": { "version": "3.0.0", "resolved": "https://registry.npmjs.org/common-path-prefix/-/common-path-prefix-3.0.0.tgz", @@ -5004,7 +4961,8 @@ "resolved": "https://registry.npmjs.org/es-module-lexer/-/es-module-lexer-1.6.0.tgz", "integrity": "sha512-qqnD1yMU6tk/jnaMosogGySTZP8YtUgAffA9nMN+E/rjxcfRQ6IEk7IiozUjgxKoFHBGjTLnrHB/YC45r/59EQ==", "dev": true, - "license": "MIT" + "license": "MIT", + "peer": true }, "node_modules/esbuild": { "version": "0.25.1", @@ -5114,88 +5072,6 @@ "url": "https://opencollective.com/eslint" } }, - "node_modules/eslint-plugin-jsdoc": { - "version": "48.11.0", - "resolved": "https://registry.npmjs.org/eslint-plugin-jsdoc/-/eslint-plugin-jsdoc-48.11.0.tgz", - "integrity": "sha512-d12JHJDPNo7IFwTOAItCeJY1hcqoIxE0lHA8infQByLilQ9xkqrRa6laWCnsuCrf+8rUnvxXY1XuTbibRBNylA==", - "dev": true, - "license": "BSD-3-Clause", - "dependencies": { - "@es-joy/jsdoccomment": "~0.46.0", - "are-docs-informative": "^0.0.2", - "comment-parser": "1.4.1", - "debug": "^4.3.5", - "escape-string-regexp": "^4.0.0", - "espree": "^10.1.0", - "esquery": "^1.6.0", - "parse-imports": "^2.1.1", - "semver": "^7.6.3", - "spdx-expression-parse": "^4.0.0", - "synckit": "^0.9.1" - }, - "engines": { - "node": ">=18" - }, - "peerDependencies": { - "eslint": "^7.0.0 || ^8.0.0 || ^9.0.0" - } - }, - "node_modules/eslint-plugin-jsdoc/node_modules/escape-string-regexp": { - "version": "4.0.0", - "resolved": "https://registry.npmjs.org/escape-string-regexp/-/escape-string-regexp-4.0.0.tgz", - "integrity": "sha512-TtpcNJ3XAzx3Gq8sWRzJaVajRs0uVxA2YAkdb1jm2YkPz4G6egUFAyA3n5vtEIZefPk5Wa4UXbKuS5fKkJWdgA==", - "dev": true, - "engines": { - "node": ">=10" - }, - "funding": { - "url": "https://github.com/sponsors/sindresorhus" - } - }, - "node_modules/eslint-plugin-jsdoc/node_modules/eslint-visitor-keys": { - "version": "4.2.0", - "resolved": "https://registry.npmjs.org/eslint-visitor-keys/-/eslint-visitor-keys-4.2.0.tgz", - "integrity": "sha512-UyLnSehNt62FFhSwjZlHmeokpRK59rcz29j+F1/aDgbkbRTk7wIc9XzdoasMUbRNKDM0qQt/+BJ4BrpFeABemw==", - "dev": true, - "license": "Apache-2.0", - "engines": { - "node": "^18.18.0 || ^20.9.0 || >=21.1.0" - }, - "funding": { - "url": "https://opencollective.com/eslint" - } - }, - "node_modules/eslint-plugin-jsdoc/node_modules/espree": { - "version": "10.3.0", - "resolved": "https://registry.npmjs.org/espree/-/espree-10.3.0.tgz", - "integrity": "sha512-0QYC8b24HWY8zjRnDTL6RiHfDbAWn63qb4LMj1Z4b076A4une81+z03Kg7l7mn/48PUTqoLptSXez8oknU8Clg==", - "dev": true, - "license": "BSD-2-Clause", - "dependencies": { - "acorn": "^8.14.0", - "acorn-jsx": "^5.3.2", - "eslint-visitor-keys": "^4.2.0" - }, - "engines": { - "node": "^18.18.0 || ^20.9.0 || >=21.1.0" - }, - "funding": { - "url": "https://opencollective.com/eslint" - } - }, - "node_modules/eslint-plugin-jsdoc/node_modules/semver": { - "version": "7.7.1", - "resolved": "https://registry.npmjs.org/semver/-/semver-7.7.1.tgz", - "integrity": "sha512-hlq8tAfn0m/61p4BVRcPzIGr6LKiMwo4VM6dGi6pt4qcRkmNzTcWq6eCEjEh+qXjkMDvPlOFFSGwQjoEa6gyMA==", - "dev": true, - "license": "ISC", - "bin": { - "semver": "bin/semver.js" - }, - "engines": { - "node": ">=10" - } - }, "node_modules/eslint-scope": { "version": "5.1.1", "resolved": "https://registry.npmjs.org/eslint-scope/-/eslint-scope-5.1.1.tgz", @@ -7859,16 +7735,6 @@ "js-yaml": "bin/js-yaml.js" } }, - "node_modules/jsdoc-type-pratt-parser": { - "version": "4.0.0", - "resolved": "https://registry.npmjs.org/jsdoc-type-pratt-parser/-/jsdoc-type-pratt-parser-4.0.0.tgz", - "integrity": "sha512-YtOli5Cmzy3q4dP26GraSOeAhqecewG04hoO8DY56CH4KJ9Fvv5qKWUCCo3HZob7esJQHCv6/+bnTy72xZZaVQ==", - "dev": true, - "license": "MIT", - "engines": { - "node": ">=12.0.0" - } - }, "node_modules/jsesc": { "version": "3.1.0", "resolved": "https://registry.npmjs.org/jsesc/-/jsesc-3.1.0.tgz", @@ -8330,20 +8196,6 @@ "node": ">= 0.10" } }, - "node_modules/parse-imports": { - "version": "2.2.1", - "resolved": "https://registry.npmjs.org/parse-imports/-/parse-imports-2.2.1.tgz", - "integrity": "sha512-OL/zLggRp8mFhKL0rNORUTR4yBYujK/uU+xZL+/0Rgm2QE4nLO9v8PzEweSJEbMGKmDRjJE4R3IMJlL2di4JeQ==", - "dev": true, - "license": "Apache-2.0 AND MIT", - "dependencies": { - "es-module-lexer": "^1.5.3", - "slashes": "^3.0.12" - }, - "engines": { - "node": ">= 18" - } - }, "node_modules/parse-json": { "version": "5.2.0", "resolved": "https://registry.npmjs.org/parse-json/-/parse-json-5.2.0.tgz", @@ -8906,13 +8758,6 @@ "node": ">=8" } }, - "node_modules/slashes": { - "version": "3.0.12", - "resolved": "https://registry.npmjs.org/slashes/-/slashes-3.0.12.tgz", - "integrity": "sha512-Q9VME8WyGkc7pJf6QEkj3wE+2CnvZMI+XJhwdTPR8Z/kWQRXi7boAWLDibRPyHRTUTPx5FaU7MsyrjI3yLB4HA==", - "dev": true, - "license": "ISC" - }, "node_modules/source-map": { "version": "0.6.1", "resolved": "https://registry.npmjs.org/source-map/-/source-map-0.6.1.tgz", @@ -8934,28 +8779,6 @@ "source-map": "^0.6.0" } }, - "node_modules/spdx-exceptions": { - "version": "2.3.0", - "resolved": "https://registry.npmjs.org/spdx-exceptions/-/spdx-exceptions-2.3.0.tgz", - "integrity": "sha512-/tTrYOC7PPI1nUAgx34hUpqXuyJG+DTHJTnIULG4rDygi4xu/tfgmq1e1cIRwRzwZgo4NLySi+ricLkZkw4i5A==", - "dev": true - }, - "node_modules/spdx-expression-parse": { - "version": "4.0.0", - "resolved": "https://registry.npmjs.org/spdx-expression-parse/-/spdx-expression-parse-4.0.0.tgz", - "integrity": "sha512-Clya5JIij/7C6bRR22+tnGXbc4VKlibKSVj2iHvVeX5iMW7s1SIQlqu699JkODJJIhh/pUu8L0/VLh8xflD+LQ==", - "dev": true, - "dependencies": { - "spdx-exceptions": "^2.1.0", - "spdx-license-ids": "^3.0.0" - } - }, - "node_modules/spdx-license-ids": { - "version": "3.0.16", - "resolved": "https://registry.npmjs.org/spdx-license-ids/-/spdx-license-ids-3.0.16.tgz", - "integrity": "sha512-eWN+LnM3GR6gPu35WxNgbGl8rmY1AEmoMDvL/QD6zYmPWgywxWqJWNdLGT+ke8dKNWrcYgYjPpG5gbTfghP8rw==", - "dev": true - }, "node_modules/sprintf-js": { "version": "1.0.3", "resolved": "https://registry.npmjs.org/sprintf-js/-/sprintf-js-1.0.3.tgz", @@ -9111,23 +8934,6 @@ "url": "https://github.com/sponsors/ljharb" } }, - "node_modules/synckit": { - "version": "0.9.2", - "resolved": "https://registry.npmjs.org/synckit/-/synckit-0.9.2.tgz", - "integrity": "sha512-vrozgXDQwYO72vHjUb/HnFbQx1exDjoKzqx23aXEg2a9VIg2TSFZ8FmeZpTjUCFMYw7mpX4BE2SFu8wI7asYsw==", - "dev": true, - "license": "MIT", - "dependencies": { - "@pkgr/core": "^0.1.0", - "tslib": "^2.6.2" - }, - "engines": { - "node": "^14.18.0 || >=16.0.0" - }, - "funding": { - "url": "https://opencollective.com/unts" - } - }, "node_modules/tapable": { "version": "2.2.1", "resolved": "https://registry.npmjs.org/tapable/-/tapable-2.2.1.tgz", @@ -9384,13 +9190,6 @@ } } }, - "node_modules/tslib": { - "version": "2.8.1", - "resolved": "https://registry.npmjs.org/tslib/-/tslib-2.8.1.tgz", - "integrity": "sha512-oJFu94HQb+KVduSUQL7wnpmqnfmLsOA/nAh6b6EH0wCEoK0/mPeXU6c3wKDV83MkOuHPRHtSXKKU99IBazS/2w==", - "dev": true, - "license": "0BSD" - }, "node_modules/tsx": { "version": "4.19.3", "resolved": "https://registry.npmjs.org/tsx/-/tsx-4.19.3.tgz", @@ -9473,6 +9272,13 @@ "node": ">=0.8.0" } }, + "node_modules/undici-types": { + "version": "7.16.0", + "resolved": "https://registry.npmjs.org/undici-types/-/undici-types-7.16.0.tgz", + "integrity": "sha512-Zz+aZWSj8LE6zoxD+xrjh4VfkIG8Ya6LvYkZqtUQGJPZjYl53ypCaUwWqo7eI0x66KBGeRo+mlBEkMSeSZ38Nw==", + "dev": true, + "license": "MIT" + }, "node_modules/unicode-canonical-property-names-ecmascript": { "version": "2.0.1", "resolved": "https://registry.npmjs.org/unicode-canonical-property-names-ecmascript/-/unicode-canonical-property-names-ecmascript-2.0.1.tgz", diff --git a/package.json b/package.json index f8a2ddea..6591b59b 100644 --- a/package.json +++ b/package.json @@ -39,6 +39,7 @@ "@babel/preset-env": "^7.26.9", "@babel/preset-typescript": "^7.26.0", "@types/jest": "^29.5.14", + "@types/node": "^24.10.1", "@typescript-eslint/eslint-plugin": "^7.18.0", "@typescript-eslint/parser": "^7.18.0", "auto-changelog": "^2.5.0", diff --git a/src/attitude/Attitude.ts b/src/attitude/Attitude.ts new file mode 100644 index 00000000..81946135 --- /dev/null +++ b/src/attitude/Attitude.ts @@ -0,0 +1,281 @@ +/** + * @author Theodore Kruczek + * @description Orbital Object ToolKit (ootk) is a collection of tools for working + * with satellites and other orbital objects. + * @license AGPL-3.0-or-later + * @copyright (c) 2025 Kruczek Labs LLC + * + * Orbital Object ToolKit is free software: you can redistribute it and/or modify it under the + * terms of the GNU Affero General Public License as published by the Free Software + * Foundation, either version 3 of the License, or (at your option) any later version. + * + * Orbital Object ToolKit is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; + * without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. + * See the GNU Affero General Public License for more details. + * + * You should have received a copy of the GNU Affero General Public License along with + * Orbital Object ToolKit. If not, see . + */ + +import { EpochUTC, Matrix, Quaternion, Vector3D } from '../main.js'; +import { Degrees, Radians, Seconds } from '../types/types.js'; +import { DEG2RAD, RAD2DEG } from '../utils/constants.js'; + +/** + * Euler angles representation (3-2-1 sequence). + */ +export interface AttitudeEulerAngles { + /** Roll angle in radians */ + roll: Radians; + /** Pitch angle in radians */ + pitch: Radians; + /** Yaw angle in radians */ + yaw: Radians; +} + +/** + * Attitude state including orientation and angular velocity. + */ +export interface AttitudeState { + /** Epoch of the state */ + epoch: EpochUTC; + /** Orientation quaternion */ + quaternion: Quaternion; + /** Angular velocity in rad/s */ + angularVelocity: Vector3D; +} + +/** + * Attitude representation and manipulation. + * Handles spacecraft orientation using quaternions and Euler angles. + */ +export class Attitude { + private quaternion_: Quaternion; + private epoch_: EpochUTC; + + constructor(quaternion: Quaternion, epoch: EpochUTC) { + this.quaternion_ = quaternion.normalize(); + this.epoch_ = epoch; + } + + /** + * Get the orientation quaternion. + * @returns Orientation quaternion + */ + get quaternion(): Quaternion { + return this.quaternion_; + } + + /** + * Get the epoch. + * @returns Epoch + */ + get epoch(): EpochUTC { + return this.epoch_; + } + + /** + * Convert quaternion to Euler angles (3-2-1 sequence). + * @returns Euler angles + */ + toEulerAngles(): AttitudeEulerAngles { + const q = this.quaternion_; + const q0 = q.w; + const q1 = q.x; + const q2 = q.y; + const q3 = q.z; + + // 3-2-1 (yaw-pitch-roll) sequence + const roll = Math.atan2(2 * (q0 * q1 + q2 * q3), 1 - 2 * (q1 * q1 + q2 * q2)) as Radians; + const pitch = Math.asin(Math.max(-1, Math.min(1, 2 * (q0 * q2 - q3 * q1)))) as Radians; + const yaw = Math.atan2(2 * (q0 * q3 + q1 * q2), 1 - 2 * (q2 * q2 + q3 * q3)) as Radians; + + return { roll, pitch, yaw }; + } + + /** + * Convert quaternion to direction cosine matrix (DCM). + * @returns Direction cosine matrix + */ + toDCM(): Matrix { + const q = this.quaternion_; + const q0 = q.w; + const q1 = q.x; + const q2 = q.y; + const q3 = q.z; + + const m11 = 1 - 2 * (q2 * q2 + q3 * q3); + const m12 = 2 * (q1 * q2 + q0 * q3); + const m13 = 2 * (q1 * q3 - q0 * q2); + + const m21 = 2 * (q1 * q2 - q0 * q3); + const m22 = 1 - 2 * (q1 * q1 + q3 * q3); + const m23 = 2 * (q2 * q3 + q0 * q1); + + const m31 = 2 * (q1 * q3 + q0 * q2); + const m32 = 2 * (q2 * q3 - q0 * q1); + const m33 = 1 - 2 * (q1 * q1 + q2 * q2); + + return new Matrix([ + [m11, m12, m13], + [m21, m22, m23], + [m31, m32, m33], + ]); + } + + /** + * Rotate a vector from body frame to inertial frame. + * @param bodyVector - Vector in body frame + * @returns Vector in inertial frame + */ + bodyToInertial(bodyVector: Vector3D): Vector3D { + return this.quaternion_.rotateVector3D(bodyVector); + } + + /** + * Rotate a vector from inertial frame to body frame. + * @param inertialVector - Vector in inertial frame + * @returns Vector in body frame + */ + inertialToBody(inertialVector: Vector3D): Vector3D { + return this.quaternion_.conjugate().rotateVector3D(inertialVector); + } + + /** + * Create attitude from Euler angles (3-2-1 sequence). + * @param roll - Roll angle in radians + * @param pitch - Pitch angle in radians + * @param yaw - Yaw angle in radians + * @param epoch - Epoch + * @returns Attitude object + */ + static fromEulerAngles(roll: Radians, pitch: Radians, yaw: Radians, epoch: EpochUTC): Attitude { + const cr = Math.cos(roll / 2); + const sr = Math.sin(roll / 2); + const cp = Math.cos(pitch / 2); + const sp = Math.sin(pitch / 2); + const cy = Math.cos(yaw / 2); + const sy = Math.sin(yaw / 2); + + const w = cr * cp * cy + sr * sp * sy; + const x = sr * cp * cy - cr * sp * sy; + const y = cr * sp * cy + sr * cp * sy; + const z = cr * cp * sy - sr * sp * cy; + + return new Attitude(new Quaternion(x, y, z, w), epoch); + } + + /** + * Create attitude from direction cosine matrix. + * @param dcm - Direction cosine matrix + * @param epoch - Epoch + * @returns Attitude object + */ + static fromDCM(dcm: Matrix, epoch: EpochUTC): Attitude { + const trace = dcm.elements[0][0] + dcm.elements[1][1] + dcm.elements[2][2]; + + let w: number, x: number, y: number, z: number; + + if (trace > 0) { + const s = Math.sqrt(trace + 1) * 2; + + w = 0.25 * s; + x = (dcm.elements[2][1] - dcm.elements[1][2]) / s; + y = (dcm.elements[0][2] - dcm.elements[2][0]) / s; + z = (dcm.elements[1][0] - dcm.elements[0][1]) / s; + } else if (dcm.elements[0][0] > dcm.elements[1][1] && dcm.elements[0][0] > dcm.elements[2][2]) { + const s = Math.sqrt(1 + dcm.elements[0][0] - dcm.elements[1][1] - dcm.elements[2][2]) * 2; + + w = (dcm.elements[2][1] - dcm.elements[1][2]) / s; + x = 0.25 * s; + y = (dcm.elements[0][1] + dcm.elements[1][0]) / s; + z = (dcm.elements[0][2] + dcm.elements[2][0]) / s; + } else if (dcm.elements[1][1] > dcm.elements[2][2]) { + const s = Math.sqrt(1 + dcm.elements[1][1] - dcm.elements[0][0] - dcm.elements[2][2]) * 2; + + w = (dcm.elements[0][2] - dcm.elements[2][0]) / s; + x = (dcm.elements[0][1] + dcm.elements[1][0]) / s; + y = 0.25 * s; + z = (dcm.elements[1][2] + dcm.elements[2][1]) / s; + } else { + const s = Math.sqrt(1 + dcm.elements[2][2] - dcm.elements[0][0] - dcm.elements[1][1]) * 2; + + w = (dcm.elements[1][0] - dcm.elements[0][1]) / s; + x = (dcm.elements[0][2] + dcm.elements[2][0]) / s; + y = (dcm.elements[1][2] + dcm.elements[2][1]) / s; + z = 0.25 * s; + } + + return new Attitude(new Quaternion(x, y, z, w), epoch); + } + + /** + * Create nadir-pointing attitude. + * Points the -Z axis toward Earth's center. + * @param position - Satellite position vector (ECI) + * @param velocity - Satellite velocity vector (ECI) + * @param epoch - Epoch + * @returns Nadir-pointing attitude + */ + static nadirPointing(position: Vector3D, velocity: Vector3D, epoch: EpochUTC): Attitude { + // Nadir pointing: -Z axis toward Earth center + const nadir = position.normalize().scale(-1); + + // Velocity direction (approximate along-track) + const alongTrack = velocity.normalize(); + + // Cross-track (orbit normal) + const crossTrack = position.cross(velocity).normalize(); + + // Recompute along-track to ensure orthogonality + const alongTrackOrth = crossTrack.cross(nadir); + + // Create DCM with body axes aligned to orbital frame + const dcm = new Matrix([ + [alongTrackOrth.x, crossTrack.x, nadir.x], + [alongTrackOrth.y, crossTrack.y, nadir.y], + [alongTrackOrth.z, crossTrack.z, nadir.z], + ]); + + return Attitude.fromDCM(dcm, epoch); + } + + /** + * Create Sun-pointing attitude. + * Points the +X axis toward the Sun. + * @param sunVector - Vector from satellite to Sun (ECI) + * @param orbitNormal - Orbit normal vector (ECI) + * @param epoch - Epoch + * @returns Sun-pointing attitude + */ + static sunPointing(sunVector: Vector3D, orbitNormal: Vector3D, epoch: EpochUTC): Attitude { + // Point +X toward Sun + const xAxis = sunVector.normalize(); + + // Keep orbit normal in Y-Z plane + const zAxis = orbitNormal.normalize(); + + // Y axis completes right-handed system + const yAxis = zAxis.cross(xAxis).normalize(); + + // Recompute Z to ensure orthogonality + const zAxisOrth = xAxis.cross(yAxis); + + const dcm = new Matrix([ + [xAxis.x, yAxis.x, zAxisOrth.x], + [xAxis.y, yAxis.y, zAxisOrth.y], + [xAxis.z, yAxis.z, zAxisOrth.z], + ]); + + return Attitude.fromDCM(dcm, epoch); + } + + /** + * Create inertially fixed attitude. + * @param epoch - Epoch + * @returns Identity attitude (no rotation) + */ + static inertial(epoch: EpochUTC): Attitude { + return new Attitude(new Quaternion(0, 0, 0, 1), epoch); + } +} diff --git a/src/attitude/AttitudePropagator.ts b/src/attitude/AttitudePropagator.ts new file mode 100644 index 00000000..cf0bfef1 --- /dev/null +++ b/src/attitude/AttitudePropagator.ts @@ -0,0 +1,233 @@ +/** + * @author Theodore Kruczek + * @description Orbital Object ToolKit (ootk) is a collection of tools for working + * with satellites and other orbital objects. + * @license AGPL-3.0-or-later + * @copyright (c) 2025 Kruczek Labs LLC + * + * Orbital Object ToolKit is free software: you can redistribute it and/or modify it under the + * terms of the GNU Affero General Public License as published by the Free Software + * Foundation, either version 3 of the License, or (at your option) any later version. + * + * Orbital Object ToolKit is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; + * without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. + * See the GNU Affero General Public License for more details. + * + * You should have received a copy of the GNU Affero General Public License along with + * Orbital Object ToolKit. If not, see . + */ + +import { EpochUTC, Quaternion, Vector3D } from '../main.js'; +import { Seconds } from '../types/types.js'; +import { Attitude, AttitudeState } from './Attitude.js'; + +/** + * Spacecraft moment of inertia tensor (kg·m²). + */ +export interface MomentOfInertia { + /** Moment of inertia about X-axis */ + Ixx: number; + /** Moment of inertia about Y-axis */ + Iyy: number; + /** Moment of inertia about Z-axis */ + Izz: number; + /** Product of inertia (default: 0 for principal axes) */ + Ixy?: number; + /** Product of inertia (default: 0 for principal axes) */ + Ixz?: number; + /** Product of inertia (default: 0 for principal axes) */ + Iyz?: number; +} + +/** + * External torque vector (N·m). + */ +export type Torque = Vector3D; + +/** + * Attitude propagator using Euler's equations of motion. + */ +export class AttitudePropagator { + private inertia_: MomentOfInertia; + private initialState_: AttitudeState; + + constructor(initialState: AttitudeState, inertia: MomentOfInertia) { + this.initialState_ = initialState; + this.inertia_ = inertia; + } + + /** + * Propagate attitude to a new epoch. + * Uses RK4 integration of Euler's equations and quaternion kinematics. + * @param targetEpoch - Target epoch + * @param externalTorque - External torque function (optional) + * @param stepSize - Integration step size in seconds (default: 1) + * @returns Attitude state at target epoch + */ + propagate( + targetEpoch: EpochUTC, + externalTorque?: (state: AttitudeState) => Torque, + stepSize: Seconds = 1 as Seconds, + ): AttitudeState { + const dt = (targetEpoch.toJulianDate() - this.initialState_.epoch.toJulianDate()) * 86400; + + if (Math.abs(dt) < 1e-6) { + return this.initialState_; + } + + let state = this.initialState_; + const steps = Math.ceil(Math.abs(dt) / stepSize); + const h = dt / steps; + + for (let i = 0; i < steps; i++) { + state = this.stepRK4_(state, h as Seconds, externalTorque); + } + + return state; + } + + /** + * Propagate torque-free motion (no external torques). + * @param targetEpoch - Target epoch + * @param stepSize - Integration step size in seconds (default: 1) + * @returns Attitude state at target epoch + */ + propagateTorqueFree(targetEpoch: EpochUTC, stepSize: Seconds = 1 as Seconds): AttitudeState { + return this.propagate(targetEpoch, undefined, stepSize); + } + + /** + * Single RK4 integration step. + */ + private stepRK4_(state: AttitudeState, h: Seconds, torqueFunc?: (state: AttitudeState) => Torque): AttitudeState { + const k1 = this.derivative_(state, torqueFunc); + const k2 = this.derivative_(this.addDerivative_(state, k1, h / 2), torqueFunc); + const k3 = this.derivative_(this.addDerivative_(state, k2, h / 2), torqueFunc); + const k4 = this.derivative_(this.addDerivative_(state, k3, h), torqueFunc); + + // Combine derivatives + const dq = new Quaternion( + (k1.dq.x + 2 * k2.dq.x + 2 * k3.dq.x + k4.dq.x) / 6, + (k1.dq.y + 2 * k2.dq.y + 2 * k3.dq.y + k4.dq.y) / 6, + (k1.dq.z + 2 * k2.dq.z + 2 * k3.dq.z + k4.dq.z) / 6, + (k1.dq.w + 2 * k2.dq.w + 2 * k3.dq.w + k4.dq.w) / 6, + ); + + const dw = new Vector3D( + (k1.dw.x + 2 * k2.dw.x + 2 * k3.dw.x + k4.dw.x) / 6, + (k1.dw.y + 2 * k2.dw.y + 2 * k3.dw.y + k4.dw.y) / 6, + (k1.dw.z + 2 * k2.dw.z + 2 * k3.dw.z + k4.dw.z) / 6, + ); + + // Update state + const newQ = new Quaternion( + state.quaternion.x + dq.x * h, + state.quaternion.y + dq.y * h, + state.quaternion.z + dq.z * h, + state.quaternion.w + dq.w * h, + ).normalize(); + + const newW = new Vector3D( + state.angularVelocity.x + dw.x * h, + state.angularVelocity.y + dw.y * h, + state.angularVelocity.z + dw.z * h, + ); + + const newEpoch = state.epoch.roll(h); + + return { + epoch: newEpoch, + quaternion: newQ, + angularVelocity: newW, + }; + } + + /** + * Calculate state derivative (dq/dt and dw/dt). + */ + private derivative_( + state: AttitudeState, + torqueFunc?: (state: AttitudeState) => Torque, + ): { dq: Quaternion; dw: Vector3D } { + const q = state.quaternion; + const w = state.angularVelocity; + + // Quaternion derivative: dq/dt = 0.5 * q * w + const wQuat = new Quaternion(w.x, w.y, w.z, 0); + const dq = q.multiply(wQuat).scale(0.5); + + // Euler's equations: I·dω/dt = τ - ω × (I·ω) + const Ixx = this.inertia_.Ixx; + const Iyy = this.inertia_.Iyy; + const Izz = this.inertia_.Izz; + + // For principal axes (products of inertia = 0) + const Iw = new Vector3D(Ixx * w.x, Iyy * w.y, Izz * w.z); + const wCrossIw = w.cross(Iw); + + // External torque + const torque = torqueFunc ? torqueFunc(state) : new Vector3D(0, 0, 0); + + // dω/dt + const dw = new Vector3D( + (torque.x - wCrossIw.x) / Ixx, + (torque.y - wCrossIw.y) / Iyy, + (torque.z - wCrossIw.z) / Izz, + ); + + return { dq, dw }; + } + + /** + * Add derivative to state for intermediate RK4 steps. + */ + private addDerivative_(state: AttitudeState, deriv: { dq: Quaternion; dw: Vector3D }, h: number): AttitudeState { + const newQ = new Quaternion( + state.quaternion.x + deriv.dq.x * h, + state.quaternion.y + deriv.dq.y * h, + state.quaternion.z + deriv.dq.z * h, + state.quaternion.w + deriv.dq.w * h, + ).normalize(); + + const newW = new Vector3D( + state.angularVelocity.x + deriv.dw.x * h, + state.angularVelocity.y + deriv.dw.y * h, + state.angularVelocity.z + deriv.dw.z * h, + ); + + return { + epoch: state.epoch.roll(h as Seconds), + quaternion: newQ, + angularVelocity: newW, + }; + } + + /** + * Create a gravity gradient torque function. + * This torque causes spacecraft to naturally align with the local vertical. + * @param orbitalRadius - Orbital radius in meters + * @returns Torque function + */ + static gravityGradientTorque(orbitalRadius: number): (state: AttitudeState) => Torque { + const mu = 3.986004418e14; // Earth's gravitational parameter (m³/s²) + const n = Math.sqrt(mu / orbitalRadius ** 3); // Orbital rate + + return (state: AttitudeState): Torque => { + // Simplified gravity gradient torque for nadir-pointing + // τ = 3n²(Izz - Iyy)sin(θ)cos(θ) (for small angles) + const attitude = new Attitude(state.quaternion, state.epoch); + const euler = attitude.toEulerAngles(); + + const Ixx = 100; // Placeholder - should use actual inertia + const Iyy = 150; + const Izz = 200; + + // Approximate torque (simplified) + const torqueX = 3 * n * n * (Izz - Iyy) * Math.sin(euler.pitch) * Math.cos(euler.pitch); + const torqueY = 3 * n * n * (Ixx - Izz) * Math.sin(euler.roll) * Math.cos(euler.roll); + const torqueZ = 3 * n * n * (Iyy - Ixx) * Math.sin(euler.yaw) * Math.cos(euler.yaw); + + return new Vector3D(torqueX, torqueY, torqueZ); + }; + } +} diff --git a/src/attitude/index.ts b/src/attitude/index.ts new file mode 100644 index 00000000..e4d881ff --- /dev/null +++ b/src/attitude/index.ts @@ -0,0 +1,21 @@ +/** + * @author Theodore Kruczek + * @description Orbital Object ToolKit (ootk) is a collection of tools for working + * with satellites and other orbital objects. + * @license AGPL-3.0-or-later + * @copyright (c) 2025 Kruczek Labs LLC + * + * Orbital Object ToolKit is free software: you can redistribute it and/or modify it under the + * terms of the GNU Affero General Public License as published by the Free Software + * Foundation, either version 3 of the License, or (at your option) any later version. + * + * Orbital Object ToolKit is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; + * without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. + * See the GNU Affero General Public License for more details. + * + * You should have received a copy of the GNU Affero General Public License along with + * Orbital Object ToolKit. If not, see . + */ + +export * from './Attitude.js'; +export * from './AttitudePropagator.js'; diff --git a/src/conjunction/ConjunctionAnalysis.ts b/src/conjunction/ConjunctionAnalysis.ts new file mode 100644 index 00000000..88188c06 --- /dev/null +++ b/src/conjunction/ConjunctionAnalysis.ts @@ -0,0 +1,255 @@ +/** + * @author Theodore Kruczek + * @description Orbital Object ToolKit (ootk) is a collection of tools for working + * with satellites and other orbital objects. + * @license AGPL-3.0-or-later + * @copyright (c) 2025 Kruczek Labs LLC + * + * Orbital Object ToolKit is free software: you can redistribute it and/or modify it under the + * terms of the GNU Affero General Public License as published by the Free Software + * Foundation, either version 3 of the License, or (at your option) any later version. + * + * Orbital Object ToolKit is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; + * without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. + * See the GNU Affero General Public License for more details. + * + * You should have received a copy of the GNU Affero General Public License along with + * Orbital Object ToolKit. If not, see . + */ + +import { EpochUTC, Matrix, Vector3D } from '../main.js'; +import { Kilometers, KilometersPerSecond, Seconds } from '../types/types.js'; +import { Satellite } from '../objects/Satellite.js'; + +/** + * Represents a close approach event between two satellites. + */ +export interface CloseApproach { + /** Time of closest approach */ + tca: EpochUTC; + /** Miss distance in kilometers */ + missDistance: Kilometers; + /** Relative velocity at TCA in km/s */ + relativeVelocity: KilometersPerSecond; + /** Position of primary satellite at TCA */ + primaryPosition: Vector3D; + /** Position of secondary satellite at TCA */ + secondaryPosition: Vector3D; + /** Velocity of primary satellite at TCA */ + primaryVelocity: Vector3D; + /** Velocity of secondary satellite at TCA */ + secondaryVelocity: Vector3D; +} + +/** + * Options for conjunction analysis. + */ +export interface ConjunctionOptions { + /** Start epoch for analysis */ + startEpoch: EpochUTC; + /** Duration to analyze in seconds */ + duration: Seconds; + /** Initial time step for coarse search in seconds (default: 60) */ + coarseStepSize?: Seconds; + /** Maximum miss distance to consider a conjunction in kilometers (default: 10) */ + maxMissDistance?: Kilometers; + /** Whether to refine TCA using interpolation (default: true) */ + refineTca?: boolean; +} + +/** + * Conjunction analyzer for detecting close approaches between satellites. + */ +export class ConjunctionAnalysis { + private primary_: Satellite; + private secondary_: Satellite; + + constructor(primary: Satellite, secondary: Satellite) { + this.primary_ = primary; + this.secondary_ = secondary; + } + + /** + * Find all close approaches during the specified time period. + * @param options - Conjunction analysis options + * @returns Array of close approach events + */ + findConjunctions(options: ConjunctionOptions): CloseApproach[] { + const { + startEpoch, + duration, + coarseStepSize = 60 as Seconds, + maxMissDistance = 10 as Kilometers, + refineTca = true, + } = options; + + const conjunctions: CloseApproach[] = []; + const steps = Math.floor(duration / coarseStepSize); + + let lastDistance: Kilometers | null = null; + let lastEpoch: EpochUTC | null = null; + let decreasingDistance = false; + + // Coarse search for local minima in distance + for (let i = 0; i <= steps; i++) { + const offset = (i * coarseStepSize) as Seconds; + const epoch = startEpoch.roll(offset); + + const distance = this.calculateDistance_(epoch); + + if (distance === null) { + lastDistance = null; + continue; + } + + if (lastDistance !== null && lastEpoch !== null) { + if (distance < lastDistance) { + decreasingDistance = true; + } else if (decreasingDistance && distance > lastDistance) { + // Found a local minimum - potential TCA + if (lastDistance <= maxMissDistance) { + let tca = lastEpoch; + let missDistance = lastDistance; + + // Refine TCA if requested + if (refineTca) { + const refined = this.refineTca_(lastEpoch, coarseStepSize); + + tca = refined.tca; + missDistance = refined.missDistance; + } + + const approach = this.createCloseApproach_(tca, missDistance); + + if (approach) { + conjunctions.push(approach); + } + } + decreasingDistance = false; + } + } + + lastDistance = distance; + lastEpoch = epoch; + } + + return conjunctions; + } + + /** + * Find the next close approach after a given epoch. + * @param startEpoch - Epoch to start searching from + * @param maxSearchDuration - Maximum time to search in seconds (default: 7 days) + * @param maxMissDistance - Maximum miss distance in kilometers (default: 10) + * @returns Next close approach or null if none found + */ + findNextConjunction( + startEpoch: EpochUTC, + maxSearchDuration: Seconds = 604800 as Seconds, + maxMissDistance: Kilometers = 10 as Kilometers, + ): CloseApproach | null { + const conjunctions = this.findConjunctions({ + startEpoch, + duration: maxSearchDuration, + maxMissDistance, + }); + + return conjunctions.length > 0 ? conjunctions[0] : null; + } + + /** + * Calculate the current distance between the two satellites. + * @param epoch - Epoch to calculate distance + * @returns Distance in kilometers or null if propagation fails + */ + private calculateDistance_(epoch: EpochUTC): Kilometers | null { + const state1 = this.primary_.toJ2000(epoch.toDateTime()); + const state2 = this.secondary_.toJ2000(epoch.toDateTime()); + + if (!state1 || !state2) { + return null; + } + + const dx = state2.position.x - state1.position.x; + const dy = state2.position.y - state1.position.y; + const dz = state2.position.z - state1.position.z; + + return Math.sqrt(dx * dx + dy * dy + dz * dz) as Kilometers; + } + + /** + * Refine the time of closest approach using golden section search. + * @param approximateTca - Approximate TCA from coarse search + * @param searchWindow - Window around approximate TCA to search (in seconds) + * @returns Refined TCA and miss distance + */ + private refineTca_( + approximateTca: EpochUTC, + searchWindow: Seconds, + ): { tca: EpochUTC; missDistance: Kilometers } { + const phi = (1 + Math.sqrt(5)) / 2; // Golden ratio + const resphi = 2 - phi; + + let a = -searchWindow / 2; + let b = searchWindow / 2; + + let x1 = a + resphi * (b - a); + let x2 = b - resphi * (b - a); + + let f1 = this.calculateDistance_(approximateTca.roll(x1 as Seconds)); + let f2 = this.calculateDistance_(approximateTca.roll(x2 as Seconds)); + + const tolerance = 0.1; // 0.1 second tolerance + + while (Math.abs(b - a) > tolerance) { + if (f1 !== null && f2 !== null && f1 < f2) { + b = x2; + x2 = x1; + f2 = f1; + x1 = a + resphi * (b - a); + f1 = this.calculateDistance_(approximateTca.roll(x1 as Seconds)); + } else { + a = x1; + x1 = x2; + f1 = f2; + x2 = b - resphi * (b - a); + f2 = this.calculateDistance_(approximateTca.roll(x2 as Seconds)); + } + } + + const tcaOffset = ((a + b) / 2) as Seconds; + const tca = approximateTca.roll(tcaOffset); + const missDistance = this.calculateDistance_(tca) ?? (0 as Kilometers); + + return { tca, missDistance }; + } + + /** + * Create a close approach object with full state information. + */ + private createCloseApproach_(tca: EpochUTC, missDistance: Kilometers): CloseApproach | null { + const state1 = this.primary_.toJ2000(tca.toDateTime()); + const state2 = this.secondary_.toJ2000(tca.toDateTime()); + + if (!state1 || !state2) { + return null; + } + + const relVel = new Vector3D( + (state2.velocity.x - state1.velocity.x) as KilometersPerSecond, + (state2.velocity.y - state1.velocity.y) as KilometersPerSecond, + (state2.velocity.z - state1.velocity.z) as KilometersPerSecond, + ); + const relSpeed = relVel.magnitude() as KilometersPerSecond; + + return { + tca, + missDistance, + relativeVelocity: relSpeed, + primaryPosition: state1.position, + secondaryPosition: state2.position, + primaryVelocity: state1.velocity, + secondaryVelocity: state2.velocity, + }; + } +} diff --git a/src/conjunction/ProbabilityOfCollision.ts b/src/conjunction/ProbabilityOfCollision.ts new file mode 100644 index 00000000..9326cab3 --- /dev/null +++ b/src/conjunction/ProbabilityOfCollision.ts @@ -0,0 +1,255 @@ +/** + * @author Theodore Kruczek + * @description Orbital Object ToolKit (ootk) is a collection of tools for working + * with satellites and other orbital objects. + * @license AGPL-3.0-or-later + * @copyright (c) 2025 Kruczek Labs LLC + * + * Orbital Object ToolKit is free software: you can redistribute it and/or modify it under the + * terms of the GNU Affero General Public License as published by the Free Software + * Foundation, either version 3 of the License, or (at your option) any later version. + * + * Orbital Object ToolKit is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; + * without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. + * See the GNU Affero General Public License for more details. + * + * You should have received a copy of the GNU Affero General Public License along with + * Orbital Object ToolKit. If not, see . + */ + +import { Matrix, StateVector, Vector3D } from '../main.js'; +import { Kilometers, KilometersPerSecond, Meters } from '../types/types.js'; + +/** + * Options for probability of collision calculation. + */ +export interface ProbabilityOfCollisionOptions { + /** Combined hard body radius (sum of both object radii) in meters */ + combinedRadius: Meters; + /** Primary object covariance matrix (3x3 position covariance in km²) */ + primaryCovariance: Matrix; + /** Secondary object covariance matrix (3x3 position covariance in km²) */ + secondaryCovariance: Matrix; + /** Method to use for calculation (default: 'chan') */ + method?: 'chan' | 'foster' | 'patera'; +} + +/** + * Result of probability of collision calculation. + */ +export interface ProbabilityOfCollisionResult { + /** Probability of collision (0-1) */ + probability: number; + /** Miss distance in kilometers */ + missDistance: Kilometers; + /** Mahalanobis distance (dimensionless) */ + mahalanobisDistance: number; + /** Combined covariance matrix */ + combinedCovariance: Matrix; + /** Method used for calculation */ + method: string; +} + +/** + * Calculates probability of collision between two space objects. + * Implements various methods including Chan, Foster, and Patera. + */ +export class ProbabilityOfCollision { + /** + * Calculate probability of collision using specified method. + * @param primaryState - State vector of primary object at TCA + * @param secondaryState - State vector of secondary object at TCA + * @param options - Calculation options + * @returns Probability of collision result + */ + static calculate( + primaryState: StateVector, + secondaryState: StateVector, + options: ProbabilityOfCollisionOptions, + ): ProbabilityOfCollisionResult { + const { combinedRadius, primaryCovariance, secondaryCovariance, method = 'chan' } = options; + + // Calculate relative state + const relPos = new Vector3D( + (secondaryState.position.x - primaryState.position.x) as Kilometers, + (secondaryState.position.y - primaryState.position.y) as Kilometers, + (secondaryState.position.z - primaryState.position.z) as Kilometers, + ); + + const relVel = new Vector3D( + (secondaryState.velocity.x - primaryState.velocity.x) as KilometersPerSecond, + (secondaryState.velocity.y - primaryState.velocity.y) as KilometersPerSecond, + (secondaryState.velocity.z - primaryState.velocity.z) as KilometersPerSecond, + ); + + const missDistance = relPos.magnitude() as Kilometers; + + // Combine covariances + const combinedCov = primaryCovariance.add(secondaryCovariance); + + // Transform to encounter frame (relative velocity aligned with x-axis) + const encounterFrame = this.createEncounterFrame_(relPos, relVel); + const covEncounter = this.transformCovariance_(combinedCov, encounterFrame); + + // Calculate probability based on method + let probability: number; + + switch (method) { + case 'foster': + probability = this.calculateFoster_(covEncounter, combinedRadius); + break; + case 'patera': + probability = this.calculatePatera_(covEncounter, combinedRadius); + break; + case 'chan': + default: + probability = this.calculateChan_(covEncounter, combinedRadius); + break; + } + + // Calculate Mahalanobis distance + const mahalanobisDistance = this.calculateMahalanobisDistance_(relPos, combinedCov); + + return { + probability, + missDistance, + mahalanobisDistance, + combinedCovariance: combinedCov, + method, + }; + } + + /** + * Calculate Pc using Chan's method (2D circular approximation). + * This is a fast, conservative approximation. + */ + private static calculateChan_(covEncounter: Matrix, combinedRadius: Meters): number { + // Project to 2D (y-z plane in encounter frame) + const sigma_y = Math.sqrt(covEncounter.elements[1][1]); + const sigma_z = Math.sqrt(covEncounter.elements[2][2]); + + // Use circular approximation + const sigma = Math.sqrt((sigma_y ** 2 + sigma_z ** 2) / 2); + const radiusKm = combinedRadius / 1000; + + // Chan's approximation + const lambda = radiusKm / sigma; + const pc = 1 - Math.exp(-(lambda ** 2) / 2); + + return Math.min(Math.max(pc, 0), 1); // Clamp to [0, 1] + } + + /** + * Calculate Pc using Foster's method (2D elliptical). + * More accurate than Chan for elliptical covariances. + */ + private static calculateFoster_(covEncounter: Matrix, combinedRadius: Meters): number { + // Extract 2D covariance in y-z plane + const Cyy = covEncounter.elements[1][1]; + const Czz = covEncounter.elements[2][2]; + const Cyz = covEncounter.elements[1][2]; + + const radiusKm = combinedRadius / 1000; + + // Calculate eigenvalues of 2D covariance + const trace = Cyy + Czz; + const det = Cyy * Czz - Cyz ** 2; + const discriminant = trace ** 2 - 4 * det; + + if (discriminant < 0 || det <= 0) { + return 0; + } + + const lambda1 = (trace + Math.sqrt(discriminant)) / 2; + const lambda2 = (trace - Math.sqrt(discriminant)) / 2; + + const sigma1 = Math.sqrt(lambda1); + const sigma2 = Math.sqrt(lambda2); + + // Foster's formula (simplified) + const alpha = radiusKm / sigma1; + const beta = sigma2 / sigma1; + + // Approximate using series expansion + const pc = 1 - Math.exp(-((alpha ** 2) / (2 * beta ** 2))); + + return Math.min(Math.max(pc, 0), 1); + } + + /** + * Calculate Pc using Patera's method (numerical integration). + * Most accurate but computationally expensive. + */ + private static calculatePatera_(covEncounter: Matrix, combinedRadius: Meters): number { + // This is a simplified version - full Patera method requires numerical integration + // For now, use Chan's method as a placeholder + // TODO: Implement full Patera numerical integration + return this.calculateChan_(covEncounter, combinedRadius); + } + + /** + * Create encounter reference frame. + * X-axis aligned with relative velocity vector. + * Y-Z plane is the collision plane. + */ + private static createEncounterFrame_(relPos: Vector3D, relVel: Vector3D): Matrix { + // X-axis: unit relative velocity + const xAxis = relVel.normalize(); + + // Z-axis: perpendicular to relative position and velocity + const zAxis = relPos.cross(relVel).normalize(); + + // Y-axis: completes right-handed system + const yAxis = zAxis.cross(xAxis); + + // Create rotation matrix (rows are the new axes in old frame) + return new Matrix([ + [xAxis.x, xAxis.y, xAxis.z], + [yAxis.x, yAxis.y, yAxis.z], + [zAxis.x, zAxis.y, zAxis.z], + ]); + } + + /** + * Transform covariance matrix to encounter frame. + */ + private static transformCovariance_(covariance: Matrix, rotation: Matrix): Matrix { + // C' = R * C * R^T + return rotation.multiply(covariance).multiply(rotation.transpose()); + } + + /** + * Calculate Mahalanobis distance. + * This is a dimensionless measure of how many standard deviations + * the miss distance is from the mean. + */ + private static calculateMahalanobisDistance_(relPos: Vector3D, covariance: Matrix): number { + const r = new Matrix([[relPos.x], [relPos.y], [relPos.z]]); + + try { + const covInv = covariance.inverse(); + const d2 = r.transpose().multiply(covInv).multiply(r).elements[0][0]; + + return Math.sqrt(Math.max(d2, 0)); + } catch { + return 0; + } + } + + /** + * Assess collision risk based on probability. + * @param probability - Probability of collision (0-1) + * @returns Risk level string + */ + static assessRisk(probability: number): string { + if (probability > 1e-4) { + return 'HIGH'; + } else if (probability > 1e-5) { + return 'MEDIUM'; + } else if (probability > 1e-6) { + return 'LOW'; + } else { + return 'NEGLIGIBLE'; + } + } +} diff --git a/src/conjunction/index.ts b/src/conjunction/index.ts new file mode 100644 index 00000000..3f3788da --- /dev/null +++ b/src/conjunction/index.ts @@ -0,0 +1,21 @@ +/** + * @author Theodore Kruczek + * @description Orbital Object ToolKit (ootk) is a collection of tools for working + * with satellites and other orbital objects. + * @license AGPL-3.0-or-later + * @copyright (c) 2025 Kruczek Labs LLC + * + * Orbital Object ToolKit is free software: you can redistribute it and/or modify it under the + * terms of the GNU Affero General Public License as published by the Free Software + * Foundation, either version 3 of the License, or (at your option) any later version. + * + * Orbital Object ToolKit is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; + * without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. + * See the GNU Affero General Public License for more details. + * + * You should have received a copy of the GNU Affero General Public License along with + * Orbital Object ToolKit. If not, see . + */ + +export * from './ConjunctionAnalysis.js'; +export * from './ProbabilityOfCollision.js'; diff --git a/src/constellation/WalkerConstellation.ts b/src/constellation/WalkerConstellation.ts new file mode 100644 index 00000000..93571bb9 --- /dev/null +++ b/src/constellation/WalkerConstellation.ts @@ -0,0 +1,304 @@ +/** + * @author Theodore Kruczek + * @description Orbital Object ToolKit (ootk) is a collection of tools for working + * with satellites and other orbital objects. + * @license AGPL-3.0-or-later + * @copyright (c) 2025 Kruczek Labs LLC + * + * Orbital Object ToolKit is free software: you can redistribute it and/or modify it under the + * terms of the GNU Affero General Public License as published by the Free Software + * Foundation, either version 3 of the License, or (at your option) any later version. + * + * Orbital Object ToolKit is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; + * without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. + * See the GNU Affero General Public License for more details. + * + * You should have received a copy of the GNU Affero General Public License along with + * Orbital Object ToolKit. If not, see . + */ + +import { ClassicalElements, EpochUTC } from '../main.js'; +import { Degrees, Kilometers, Radians } from '../types/types.js'; +import { DEG2RAD, RAD2DEG } from '../utils/constants.js'; +import { Satellite } from '../objects/Satellite.js'; + +/** + * Walker constellation parameters. + * Walker constellations are denoted as i:t/p/f where: + * - i = inclination + * - t = total number of satellites + * - p = number of equally spaced orbital planes + * - f = relative spacing between satellites in adjacent planes (phasing parameter) + */ +export interface WalkerConstellationParams { + /** Total number of satellites */ + totalSatellites: number; + /** Number of orbital planes */ + numberOfPlanes: number; + /** Relative phasing between planes (0 to numberOfPlanes-1) */ + relativePhasing: number; + /** Orbital inclination in degrees */ + inclination: Degrees; + /** Semi-major axis in kilometers */ + semiMajorAxis: Kilometers; + /** Eccentricity (default: 0 for circular) */ + eccentricity?: number; + /** Argument of periapsis in degrees (default: 0) */ + argumentOfPeriapsis?: Degrees; + /** Epoch for the constellation */ + epoch: EpochUTC; +} + +/** + * Represents a satellite in a Walker constellation. + */ +export interface WalkerSatellite { + /** Satellite index (0 to totalSatellites-1) */ + index: number; + /** Plane number (0 to numberOfPlanes-1) */ + plane: number; + /** Position in plane (0 to satellitesPerPlane-1) */ + positionInPlane: number; + /** Classical orbital elements */ + elements: ClassicalElements; + /** Satellite object (if TLE creation is enabled) */ + satellite?: Satellite; +} + +/** + * Walker constellation generator. + * Creates satellite constellations using the Walker Delta pattern. + */ +export class WalkerConstellation { + private params_: WalkerConstellationParams; + private satellites_: WalkerSatellite[] = []; + + constructor(params: WalkerConstellationParams) { + this.validateParams_(params); + this.params_ = params; + this.generateConstellation_(); + } + + /** + * Get all satellites in the constellation. + * @returns Array of Walker satellites + */ + getSatellites(): WalkerSatellite[] { + return this.satellites_; + } + + /** + * Get satellites in a specific orbital plane. + * @param planeNumber - Plane number (0 to numberOfPlanes-1) + * @returns Array of satellites in the plane + */ + getSatellitesInPlane(planeNumber: number): WalkerSatellite[] { + return this.satellites_.filter((sat) => sat.plane === planeNumber); + } + + /** + * Get the orbital elements for a specific satellite. + * @param index - Satellite index (0 to totalSatellites-1) + * @returns Classical orbital elements or undefined if index invalid + */ + getSatelliteElements(index: number): ClassicalElements | undefined { + return this.satellites_[index]?.elements; + } + + /** + * Get the Walker constellation designation string. + * @returns Walker designation (e.g., "53:6/6/1") + */ + getDesignation(): string { + const { inclination, totalSatellites, numberOfPlanes, relativePhasing } = this.params_; + + return `${inclination}:${totalSatellites}/${numberOfPlanes}/${relativePhasing}`; + } + + /** + * Calculate the street-of-coverage width for the constellation. + * @param minElevation - Minimum elevation angle in degrees + * @returns Street width in degrees + */ + calculateStreetWidth(minElevation: Degrees): Degrees { + const { semiMajorAxis } = this.params_; + const Re = 6371; // Earth radius in km + const h = semiMajorAxis - Re; + + // Calculate nadir angle + const elRad = minElevation * DEG2RAD; + const rho = Math.acos((Re / (Re + h)) * Math.cos(elRad)); + const eta = (Math.PI / 2 - elRad - rho) as Radians; + + // Angular radius of coverage + const lambda = (eta * RAD2DEG) as Degrees; + + return (2 * lambda) as Degrees; + } + + /** + * Validate Walker constellation parameters. + */ + private validateParams_(params: WalkerConstellationParams): void { + if (params.totalSatellites <= 0) { + throw new Error('Total satellites must be positive'); + } + + if (params.numberOfPlanes <= 0) { + throw new Error('Number of planes must be positive'); + } + + if (params.totalSatellites % params.numberOfPlanes !== 0) { + throw new Error('Total satellites must be divisible by number of planes'); + } + + if (params.relativePhasing < 0 || params.relativePhasing >= params.numberOfPlanes) { + throw new Error('Relative phasing must be between 0 and numberOfPlanes-1'); + } + + if (params.inclination < 0 || params.inclination > 180) { + throw new Error('Inclination must be between 0 and 180 degrees'); + } + + if (params.semiMajorAxis <= 6371) { + throw new Error('Semi-major axis must be greater than Earth radius (6371 km)'); + } + } + + /** + * Generate the constellation satellites. + */ + private generateConstellation_(): void { + const { + totalSatellites, + numberOfPlanes, + relativePhasing, + inclination, + semiMajorAxis, + eccentricity = 0, + argumentOfPeriapsis = 0 as Degrees, + epoch, + } = this.params_; + + const satellitesPerPlane = totalSatellites / numberOfPlanes; + const raanSpacing = (360 / numberOfPlanes) as Degrees; + const taSpacing = (360 / satellitesPerPlane) as Degrees; + + let satIndex = 0; + + for (let plane = 0; plane < numberOfPlanes; plane++) { + const raan = (plane * raanSpacing) as Degrees; + + for (let pos = 0; pos < satellitesPerPlane; pos++) { + // Calculate true anomaly with phasing + const phaseOffset = (plane * relativePhasing * taSpacing) / numberOfPlanes; + const ta = ((pos * taSpacing + phaseOffset) % 360) as Degrees; + + const elements = new ClassicalElements({ + epoch, + semimajorAxis: semiMajorAxis, + eccentricity, + inclination: (inclination * DEG2RAD) as Radians, + rightAscension: (raan * DEG2RAD) as Radians, + argPerigee: (argumentOfPeriapsis * DEG2RAD) as Radians, + trueAnomaly: (ta * DEG2RAD) as Radians, + }); + + this.satellites_.push({ + index: satIndex, + plane, + positionInPlane: pos, + elements, + }); + + satIndex++; + } + } + } + + /** + * Create a Walker Delta constellation (most common pattern). + * @param t - Total number of satellites + * @param p - Number of orbital planes + * @param f - Relative phasing parameter + * @param inclination - Orbital inclination in degrees + * @param altitude - Orbital altitude in kilometers + * @param epoch - Epoch for the constellation + * @returns WalkerConstellation instance + */ + static createDelta( + t: number, + p: number, + f: number, + inclination: Degrees, + altitude: Kilometers, + epoch: EpochUTC, + ): WalkerConstellation { + const Re = 6371; // Earth radius + const semiMajorAxis = (Re + altitude) as Kilometers; + + return new WalkerConstellation({ + totalSatellites: t, + numberOfPlanes: p, + relativePhasing: f, + inclination, + semiMajorAxis, + eccentricity: 0, + epoch, + }); + } + + /** + * Create a Walker Star constellation (polar, counter-rotating planes). + * @param t - Total number of satellites + * @param p - Number of orbital planes + * @param f - Relative phasing parameter + * @param altitude - Orbital altitude in kilometers + * @param epoch - Epoch for the constellation + * @returns WalkerConstellation instance + */ + static createStar(t: number, p: number, f: number, altitude: Kilometers, epoch: EpochUTC): WalkerConstellation { + // Star pattern uses polar orbits (90 degrees) + return WalkerConstellation.createDelta(t, p, f, 90 as Degrees, altitude, epoch); + } + + /** + * Create a preset GPS-like constellation. + * GPS uses 24 satellites in 6 planes at 55° inclination. + * @param epoch - Epoch for the constellation + * @returns WalkerConstellation instance + */ + static createGPS(epoch: EpochUTC): WalkerConstellation { + return WalkerConstellation.createDelta(24, 6, 1, 55 as Degrees, 20180 as Kilometers, epoch); + } + + /** + * Create a preset Galileo-like constellation. + * Galileo uses 24 satellites in 3 planes at 56° inclination. + * @param epoch - Epoch for the constellation + * @returns WalkerConstellation instance + */ + static createGalileo(epoch: EpochUTC): WalkerConstellation { + return WalkerConstellation.createDelta(24, 3, 1, 56 as Degrees, 23222 as Kilometers, epoch); + } + + /** + * Create a preset Iridium-like constellation. + * Iridium uses 66 satellites in 6 planes at 86.4° inclination. + * @param epoch - Epoch for the constellation + * @returns WalkerConstellation instance + */ + static createIridium(epoch: EpochUTC): WalkerConstellation { + return WalkerConstellation.createDelta(66, 6, 1, 86.4 as Degrees, 780 as Kilometers, epoch); + } + + /** + * Create a preset OneWeb-like constellation. + * OneWeb uses 648 satellites in 18 planes at 87.9° inclination. + * @param epoch - Epoch for the constellation + * @returns WalkerConstellation instance + */ + static createOneWeb(epoch: EpochUTC): WalkerConstellation { + return WalkerConstellation.createDelta(648, 18, 1, 87.9 as Degrees, 1200 as Kilometers, epoch); + } +} diff --git a/src/constellation/index.ts b/src/constellation/index.ts new file mode 100644 index 00000000..c7b01cd0 --- /dev/null +++ b/src/constellation/index.ts @@ -0,0 +1,20 @@ +/** + * @author Theodore Kruczek + * @description Orbital Object ToolKit (ootk) is a collection of tools for working + * with satellites and other orbital objects. + * @license AGPL-3.0-or-later + * @copyright (c) 2025 Kruczek Labs LLC + * + * Orbital Object ToolKit is free software: you can redistribute it and/or modify it under the + * terms of the GNU Affero General Public License as published by the Free Software + * Foundation, either version 3 of the License, or (at your option) any later version. + * + * Orbital Object ToolKit is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; + * without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. + * See the GNU Affero General Public License for more details. + * + * You should have received a copy of the GNU Affero General Public License along with + * Orbital Object ToolKit. If not, see . + */ + +export * from './WalkerConstellation.js'; diff --git a/src/coverage/AccessStatistics.ts b/src/coverage/AccessStatistics.ts new file mode 100644 index 00000000..5d8614d4 --- /dev/null +++ b/src/coverage/AccessStatistics.ts @@ -0,0 +1,231 @@ +/** + * @author Theodore Kruczek + * @description Orbital Object ToolKit (ootk) is a collection of tools for working + * with satellites and other orbital objects. + * @license AGPL-3.0-or-later + * @copyright (c) 2025 Kruczek Labs LLC + * + * Orbital Object ToolKit is free software: you can redistribute it and/or modify it under the + * terms of the GNU Affero General Public License as published by the Free Software + * Foundation, either version 3 of the License, or (at your option) any later version. + * + * Orbital Object ToolKit is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; + * without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. + * See the GNU Affero General Public License for more details. + * + * You should have received a copy of the GNU Affero General Public License along with + * Orbital Object ToolKit. If not, see . + */ + +import { EpochUTC } from '../main.js'; +import { Degrees, Seconds } from '../types/types.js'; +import { Sensor } from '../objects/Sensor.js'; +import { Satellite } from '../objects/Satellite.js'; + +/** + * Represents an access interval when a satellite is visible to a sensor. + */ +export interface AccessInterval { + /** Start epoch of access */ + start: EpochUTC; + /** End epoch of access */ + end: EpochUTC; + /** Duration in seconds */ + duration: Seconds; + /** Maximum elevation during this access */ + maxElevation: Degrees; + /** Time of maximum elevation */ + timeOfMaxElevation: EpochUTC; +} + +/** + * Statistics about satellite access over a time period. + */ +export interface AccessStats { + /** Total number of access intervals */ + totalAccesses: number; + /** Total access time in seconds */ + totalAccessTime: Seconds; + /** Mean access time in seconds */ + meanAccessTime: Seconds; + /** Minimum access time in seconds */ + minAccessTime: Seconds; + /** Maximum access time in seconds */ + maxAccessTime: Seconds; + /** Mean time between accesses in seconds */ + meanGapTime: Seconds; + /** Maximum time between accesses in seconds */ + maxGapTime: Seconds; + /** Coverage percentage (0-100) */ + coveragePercentage: number; + /** Mean maximum elevation in degrees */ + meanMaxElevation: Degrees; + /** All access intervals */ + intervals: AccessInterval[]; +} + +/** + * Calculates access statistics for satellite-to-ground visibility. + */ +export class AccessStatistics { + private sensor_: Sensor; + private satellite_: Satellite; + + constructor(sensor: Sensor, satellite: Satellite) { + this.sensor_ = sensor; + this.satellite_ = satellite; + } + + /** + * Calculate access statistics over a time period. + * @param startEpoch - Start of analysis period + * @param duration - Duration of analysis in seconds + * @param stepSize - Time step for analysis in seconds (default: 10) + * @returns Access statistics + */ + calculate(startEpoch: EpochUTC, duration: Seconds, stepSize: Seconds = 10 as Seconds): AccessStats { + const intervals: AccessInterval[] = []; + const steps = Math.floor(duration / stepSize); + + let inAccess = false; + let accessStart: EpochUTC | null = null; + let maxEl: Degrees = 0 as Degrees; + let timeOfMaxEl: EpochUTC | null = null; + + // Scan through time to find access intervals + for (let i = 0; i <= steps; i++) { + const offset = (i * stepSize) as Seconds; + const epoch = startEpoch.roll(offset); + + const rae = this.sensor_.rae(this.satellite_, epoch.toDateTime()); + + if (!rae) { + continue; + } + + const isVisible = this.sensor_.isRaeInFov(rae); + + if (isVisible && !inAccess) { + // Start of new access + inAccess = true; + accessStart = epoch; + maxEl = rae.el; + timeOfMaxEl = epoch; + } else if (isVisible && inAccess) { + // During access, track max elevation + if (rae.el > maxEl) { + maxEl = rae.el; + timeOfMaxEl = epoch; + } + } else if (!isVisible && inAccess) { + // End of access + inAccess = false; + if (accessStart && timeOfMaxEl) { + const accessDuration = ((epoch.toJulianDate() - accessStart.toJulianDate()) * 86400) as Seconds; + + intervals.push({ + start: accessStart, + end: epoch, + duration: accessDuration, + maxElevation: maxEl, + timeOfMaxElevation: timeOfMaxEl, + }); + } + maxEl = 0 as Degrees; + timeOfMaxEl = null; + } + } + + // Handle case where access is still ongoing at end of period + if (inAccess && accessStart && timeOfMaxEl) { + const endEpoch = startEpoch.roll(duration); + const accessDuration = ((endEpoch.toJulianDate() - accessStart.toJulianDate()) * 86400) as Seconds; + + intervals.push({ + start: accessStart, + end: endEpoch, + duration: accessDuration, + maxElevation: maxEl, + timeOfMaxElevation: timeOfMaxEl, + }); + } + + return this.computeStats(intervals, duration); + } + + /** + * Compute statistics from access intervals. + * @param intervals - Array of access intervals + * @param totalDuration - Total duration of analysis period + * @returns Computed statistics + */ + private computeStats(intervals: AccessInterval[], totalDuration: Seconds): AccessStats { + if (intervals.length === 0) { + return { + totalAccesses: 0, + totalAccessTime: 0 as Seconds, + meanAccessTime: 0 as Seconds, + minAccessTime: 0 as Seconds, + maxAccessTime: 0 as Seconds, + meanGapTime: 0 as Seconds, + maxGapTime: 0 as Seconds, + coveragePercentage: 0, + meanMaxElevation: 0 as Degrees, + intervals: [], + }; + } + + const totalAccessTime = intervals.reduce((sum, interval) => sum + interval.duration, 0) as Seconds; + const durations = intervals.map((interval) => interval.duration); + const elevations = intervals.map((interval) => interval.maxElevation); + + const minAccessTime = Math.min(...durations) as Seconds; + const maxAccessTime = Math.max(...durations) as Seconds; + const meanAccessTime = (totalAccessTime / intervals.length) as Seconds; + const meanMaxElevation = (elevations.reduce((sum, el) => sum + el, 0) / intervals.length) as Degrees; + + // Calculate gap times + const gapTimes: number[] = []; + + for (let i = 0; i < intervals.length - 1; i++) { + const gapTime = (intervals[i + 1].start.toJulianDate() - intervals[i].end.toJulianDate()) * 86400; + + gapTimes.push(gapTime); + } + + const meanGapTime = gapTimes.length > 0 ? ((gapTimes.reduce((sum, gap) => sum + gap, 0) / gapTimes.length) as Seconds) : (0 as Seconds); + const maxGapTime = gapTimes.length > 0 ? (Math.max(...gapTimes) as Seconds) : (0 as Seconds); + + const coveragePercentage = (totalAccessTime / totalDuration) * 100; + + return { + totalAccesses: intervals.length, + totalAccessTime, + meanAccessTime, + minAccessTime, + maxAccessTime, + meanGapTime, + maxGapTime, + coveragePercentage, + meanMaxElevation, + intervals, + }; + } + + /** + * Find the next access opportunity after a given epoch. + * @param startEpoch - Epoch to start searching from + * @param maxSearchDuration - Maximum time to search in seconds (default: 7 days) + * @param stepSize - Time step for search in seconds (default: 60) + * @returns Next access interval or null if none found + */ + findNextAccess( + startEpoch: EpochUTC, + maxSearchDuration: Seconds = 604800 as Seconds, + stepSize: Seconds = 60 as Seconds, + ): AccessInterval | null { + const stats = this.calculate(startEpoch, maxSearchDuration, stepSize); + + return stats.intervals.length > 0 ? stats.intervals[0] : null; + } +} diff --git a/src/coverage/AreaCoverage.ts b/src/coverage/AreaCoverage.ts new file mode 100644 index 00000000..4a1d3952 --- /dev/null +++ b/src/coverage/AreaCoverage.ts @@ -0,0 +1,279 @@ +/** + * @author Theodore Kruczek + * @description Orbital Object ToolKit (ootk) is a collection of tools for working + * with satellites and other orbital objects. + * @license AGPL-3.0-or-later + * @copyright (c) 2025 Kruczek Labs LLC + * + * Orbital Object ToolKit is free software: you can redistribute it and/or modify it under the + * terms of the GNU Affero General Public License as published by the Free Software + * Foundation, either version 3 of the License, or (at your option) any later version. + * + * Orbital Object ToolKit is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; + * without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. + * See the GNU Affero General Public License for more details. + * + * You should have received a copy of the GNU Affero General Public License along with + * Orbital Object ToolKit. If not, see . + */ + +import { EpochUTC } from '../main.js'; +import { Degrees, Kilometers, Seconds } from '../types/types.js'; +import { Satellite } from '../objects/Satellite.js'; +import { GroundObject } from '../objects/GroundObject.js'; +import { DEG2RAD, RAD2DEG, RADIUS_OF_EARTH } from '../utils/constants.js'; + +/** + * Represents a grid point for coverage analysis. + */ +export interface CoverageGridPoint { + /** Latitude in degrees */ + lat: Degrees; + /** Longitude in degrees */ + lon: Degrees; + /** Number of times this point was covered */ + coverageCount: number; + /** Total coverage time in seconds */ + coverageTime: Seconds; + /** Maximum elevation angle observed from this point */ + maxElevation: Degrees; + /** Percentage of time covered (0-100) */ + coveragePercentage: number; +} + +/** + * Options for area coverage calculation. + */ +export interface AreaCoverageOptions { + /** Start epoch */ + startEpoch: EpochUTC; + /** Duration in seconds */ + duration: Seconds; + /** Time step in seconds (default: 60) */ + stepSize?: Seconds; + /** Latitude resolution in degrees (default: 5) */ + latResolution?: Degrees; + /** Longitude resolution in degrees (default: 5) */ + lonResolution?: Degrees; + /** Minimum elevation angle in degrees (default: 0) */ + minElevation?: Degrees; + /** Latitude bounds (default: -90 to 90) */ + latBounds?: { min: Degrees; max: Degrees }; + /** Longitude bounds (default: -180 to 180) */ + lonBounds?: { min: Degrees; max: Degrees }; +} + +/** + * Area coverage analyzer for satellite coverage analysis. + * Computes which ground locations are covered by a satellite over time. + */ +export class AreaCoverage { + private satellites_: Satellite[]; + + constructor(satellites: Satellite | Satellite[]) { + this.satellites_ = Array.isArray(satellites) ? satellites : [satellites]; + } + + /** + * Calculate area coverage grid. + * @param options - Coverage calculation options + * @returns Array of coverage grid points + */ + calculateGrid(options: AreaCoverageOptions): CoverageGridPoint[] { + const { + startEpoch, + duration, + stepSize = 60 as Seconds, + latResolution = 5 as Degrees, + lonResolution = 5 as Degrees, + minElevation = 0 as Degrees, + latBounds = { min: -90 as Degrees, max: 90 as Degrees }, + lonBounds = { min: -180 as Degrees, max: 180 as Degrees }, + } = options; + + // Create coverage grid + const grid = this.initializeGrid_(latBounds, lonBounds, latResolution, lonResolution); + + const steps = Math.floor(duration / stepSize); + + // Scan through time and update coverage for each point + for (let i = 0; i <= steps; i++) { + const offset = (i * stepSize) as Seconds; + const epoch = startEpoch.roll(offset); + + for (const satellite of this.satellites_) { + const state = satellite.toJ2000(epoch.toDateTime()); + + if (!state) { + continue; + } + + const satPos = state.toITRF().toGeodetic(); + + // For each grid point, check if satellite is visible + for (const point of grid) { + const groundPos = new GroundObject({ lat: point.lat, lon: point.lon, alt: 0 as Kilometers }); + const el = this.calculateElevation_(satPos.latDeg as Degrees, satPos.lonDeg as Degrees, satPos.alt, groundPos); + + if (el >= minElevation) { + point.coverageCount++; + point.coverageTime = (point.coverageTime + stepSize) as Seconds; + point.maxElevation = Math.max(point.maxElevation, el) as Degrees; + } + } + } + } + + // Calculate coverage percentages + for (const point of grid) { + point.coveragePercentage = (point.coverageTime / duration) * 100; + } + + return grid; + } + + /** + * Calculate instantaneous coverage area at a specific time. + * @param epoch - Epoch to calculate coverage + * @param minElevation - Minimum elevation angle in degrees + * @param resolution - Grid resolution in degrees + * @returns Array of grid points that are currently covered + */ + calculateInstantaneousCoverage( + epoch: EpochUTC, + minElevation: Degrees = 0 as Degrees, + resolution: Degrees = 5 as Degrees, + ): CoverageGridPoint[] { + return this.calculateGrid({ + startEpoch: epoch, + duration: 1 as Seconds, + stepSize: 1 as Seconds, + latResolution: resolution, + lonResolution: resolution, + minElevation, + }).filter((point) => point.coverageCount > 0); + } + + /** + * Calculate coverage percentage over a region. + * @param options - Coverage calculation options + * @returns Coverage percentage (0-100) + */ + calculateCoveragePercentage(options: AreaCoverageOptions): number { + const grid = this.calculateGrid(options); + const coveredPoints = grid.filter((point) => point.coverageCount > 0).length; + + return (coveredPoints / grid.length) * 100; + } + + /** + * Calculate the instantaneous swath width at a given elevation angle. + * @param altitude - Satellite altitude in kilometers + * @param minElevation - Minimum elevation angle in degrees + * @returns Swath width in kilometers + */ + static calculateSwathWidth(altitude: Kilometers, minElevation: Degrees = 0 as Degrees): Kilometers { + const Re = RADIUS_OF_EARTH; + const h = altitude; + const elRad = minElevation * DEG2RAD; + + // Calculate nadir angle using spherical geometry + const rho = Math.acos((Re / (Re + h)) * Math.cos(elRad)); + const eta = Math.PI / 2 - elRad - rho; + + // Calculate swath half-angle + const lambda = Re * eta; + + // Total swath width + return (2 * lambda) as Kilometers; + } + + /** + * Calculate the footprint radius for a satellite. + * @param altitude - Satellite altitude in kilometers + * @param minElevation - Minimum elevation angle in degrees + * @returns Footprint radius in kilometers + */ + static calculateFootprintRadius(altitude: Kilometers, minElevation: Degrees = 0 as Degrees): Kilometers { + const swathWidth = AreaCoverage.calculateSwathWidth(altitude, minElevation); + + return (swathWidth / 2) as Kilometers; + } + + /** + * Calculate maximum coverage latitude for an inclined orbit. + * @param inclination - Orbital inclination in degrees + * @param altitude - Satellite altitude in kilometers + * @param minElevation - Minimum elevation angle in degrees + * @returns Maximum latitude that can be covered + */ + static calculateMaxCoverageLatitude( + inclination: Degrees, + altitude: Kilometers, + minElevation: Degrees = 0 as Degrees, + ): Degrees { + const footprintRadius = AreaCoverage.calculateFootprintRadius(altitude, minElevation); + const angularRadius = (footprintRadius / RADIUS_OF_EARTH) * RAD2DEG; + + return (inclination + angularRadius) as Degrees; + } + + /** + * Initialize coverage grid. + */ + private initializeGrid_( + latBounds: { min: Degrees; max: Degrees }, + lonBounds: { min: Degrees; max: Degrees }, + latResolution: Degrees, + lonResolution: Degrees, + ): CoverageGridPoint[] { + const grid: CoverageGridPoint[] = []; + + for (let lat = latBounds.min as number; lat <= latBounds.max; lat += latResolution) { + for (let lon = lonBounds.min as number; lon < lonBounds.max; lon += lonResolution) { + grid.push({ + lat: lat as Degrees, + lon: lon as Degrees, + coverageCount: 0, + coverageTime: 0 as Seconds, + maxElevation: 0 as Degrees, + coveragePercentage: 0, + }); + } + } + + return grid; + } + + /** + * Calculate elevation angle from a ground point to a satellite. + */ + private calculateElevation_( + satLat: Degrees, + satLon: Degrees, + satAlt: Kilometers, + groundPos: GroundObject, + ): Degrees { + // Convert to radians + const lat1 = groundPos.lat * DEG2RAD; + const lon1 = groundPos.lon * DEG2RAD; + const lat2 = satLat * DEG2RAD; + const lon2 = satLon * DEG2RAD; + + // Calculate angular distance + const dLat = lat2 - lat1; + const dLon = lon2 - lon1; + const a = Math.sin(dLat / 2) ** 2 + Math.cos(lat1) * Math.cos(lat2) * Math.sin(dLon / 2) ** 2; + const c = 2 * Math.atan2(Math.sqrt(a), Math.sqrt(1 - a)); + const distance = RADIUS_OF_EARTH * c; + + // Calculate elevation angle + const Re = RADIUS_OF_EARTH; + const h = satAlt; + const d = distance; + + const elevation = Math.atan2(h - Re * (1 - Math.cos(d / Re)), Re * Math.sin(d / Re)); + + return (elevation * RAD2DEG) as Degrees; + } +} diff --git a/src/coverage/GroundTrack.ts b/src/coverage/GroundTrack.ts new file mode 100644 index 00000000..5de930bb --- /dev/null +++ b/src/coverage/GroundTrack.ts @@ -0,0 +1,210 @@ +/** + * @author Theodore Kruczek + * @description Orbital Object ToolKit (ootk) is a collection of tools for working + * with satellites and other orbital objects. + * @license AGPL-3.0-or-later + * @copyright (c) 2025 Kruczek Labs LLC + * + * Orbital Object ToolKit is free software: you can redistribute it and/or modify it under the + * terms of the GNU Affero General Public License as published by the Free Software + * Foundation, either version 3 of the License, or (at your option) any later version. + * + * Orbital Object ToolKit is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; + * without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. + * See the GNU Affero General Public License for more details. + * + * You should have received a copy of the GNU Affero General Public License along with + * Orbital Object ToolKit. If not, see . + */ + +import { EpochUTC } from '../main.js'; +import { GroundObject } from '../objects/GroundObject.js'; +import { Degrees, Kilometers, Seconds } from '../types/types.js'; +import { Satellite } from '../objects/Satellite.js'; + +/** + * Represents a point on a ground track. + */ +export interface GroundTrackPoint { + /** Epoch of the point */ + epoch: EpochUTC; + /** Latitude in degrees */ + lat: Degrees; + /** Longitude in degrees */ + lon: Degrees; + /** Altitude in kilometers */ + alt: Kilometers; +} + +/** + * Options for ground track calculation. + */ +export interface GroundTrackOptions { + /** Start epoch */ + startEpoch: EpochUTC; + /** Duration in seconds */ + duration: Seconds; + /** Time step in seconds (default: 60) */ + stepSize?: Seconds; + /** Whether to include ascending/descending node crossings */ + includeNodeCrossings?: boolean; +} + +/** + * Ground track calculator for satellite orbits. + * Generates ground track points showing where a satellite passes over the Earth's surface. + */ +export class GroundTrack { + private satellite_: Satellite; + + constructor(satellite: Satellite) { + this.satellite_ = satellite; + } + + /** + * Calculate ground track points for the satellite. + * @param options - Ground track calculation options + * @returns Array of ground track points + */ + calculate(options: GroundTrackOptions): GroundTrackPoint[] { + const { startEpoch, duration, stepSize = 60 as Seconds, includeNodeCrossings = false } = options; + const points: GroundTrackPoint[] = []; + const steps = Math.floor(duration / stepSize); + + let lastLat: Degrees | null = null; + + for (let i = 0; i <= steps; i++) { + const offset = (i * stepSize) as Seconds; + const epoch = startEpoch.roll(offset); + const state = this.satellite_.toJ2000(epoch.toDateTime()); + + if (!state) { + continue; + } + + const lla = state.toITRF().toGeodetic(); + + const point: GroundTrackPoint = { + epoch, + lat: lla.latDeg as Degrees, + lon: lla.lonDeg as Degrees, + alt: lla.alt, + }; + + points.push(point); + + // Detect node crossings (equator crossings) + if (includeNodeCrossings && lastLat !== null) { + if ((lastLat < 0 && lla.latDeg >= 0) || (lastLat > 0 && lla.latDeg <= 0)) { + // Mark this as a node crossing + (point as GroundTrackPoint & { isNodeCrossing?: boolean }).isNodeCrossing = true; + } + } + + lastLat = lla.latDeg as Degrees; + } + + return points; + } + + /** + * Calculate ground track for one complete orbit. + * @param startEpoch - Starting epoch + * @param stepSize - Time step in seconds (default: 60) + * @returns Array of ground track points for one orbit + */ + calculateOrbit(startEpoch: EpochUTC, stepSize: Seconds = 60 as Seconds): GroundTrackPoint[] { + // Get orbital period + const state = this.satellite_.toJ2000(startEpoch.toDateTime()); + + if (!state) { + return []; + } + + const elements = state.toClassicalElements(); + const period = elements.period; + + return this.calculate({ + startEpoch, + duration: (period * 60) as Seconds, + stepSize, + includeNodeCrossings: true, + }); + } + + /** + * Get the subsatellite point (ground position directly below satellite). + * @param epoch - Epoch to calculate subsatellite point + * @returns Ground position or null if propagation fails + */ + getSubsatellitePoint(epoch: EpochUTC): GroundObject | null { + const state = this.satellite_.toJ2000(epoch.toDateTime()); + + if (!state) { + return null; + } + + const lla = state.toITRF().toGeodetic(); + + return new GroundObject({ lat: lla.latDeg as Degrees, lon: lla.lonDeg as Degrees, alt: lla.alt }); + } + + /** + * Calculate maximum and minimum latitudes reached by the satellite. + * For circular orbits, this equals the inclination. + * @param startEpoch - Starting epoch + * @returns Object with max and min latitudes + */ + getLatitudeBounds(startEpoch: EpochUTC): { maxLat: Degrees; minLat: Degrees } { + const state = this.satellite_.toJ2000(startEpoch.toDateTime()); + + if (!state) { + return { maxLat: 0 as Degrees, minLat: 0 as Degrees }; + } + + const elements = state.toClassicalElements(); + const inc = elements.inclinationDegrees; + + return { + maxLat: inc, + minLat: (-inc) as Degrees, + }; + } + + /** + * Calculate the repeat ground track parameters. + * A repeat ground track occurs when the satellite returns to the same ground track + * after a whole number of orbits and Earth rotations. + * @param startEpoch - Starting epoch + * @returns Repeat parameters or null if satellite doesn't have a repeat ground track + */ + getRepeatGroundTrack(startEpoch: EpochUTC): { + orbits: number; + days: number; + error: number; + } | null { + const state = this.satellite_.toJ2000(startEpoch.toDateTime()); + + if (!state) { + return null; + } + + const elements = state.toClassicalElements(); + const period = elements.period; // in seconds + const orbitsPerDay = 86400 / period; + + // Try to find a simple ratio (within 1%) + for (let days = 1; days <= 30; days++) { + const orbits = Math.round(orbitsPerDay * days); + const actualRatio = orbits / days; + const error = Math.abs(actualRatio - orbitsPerDay) / orbitsPerDay; + + if (error < 0.01) { + // Within 1% + return { orbits, days, error }; + } + } + + return null; + } +} diff --git a/src/coverage/index.ts b/src/coverage/index.ts new file mode 100644 index 00000000..d938a5e8 --- /dev/null +++ b/src/coverage/index.ts @@ -0,0 +1,22 @@ +/** + * @author Theodore Kruczek + * @description Orbital Object ToolKit (ootk) is a collection of tools for working + * with satellites and other orbital objects. + * @license AGPL-3.0-or-later + * @copyright (c) 2025 Kruczek Labs LLC + * + * Orbital Object ToolKit is free software: you can redistribute it and/or modify it under the + * terms of the GNU Affero General Public License as published by the Free Software + * Foundation, either version 3 of the License, or (at your option) any later version. + * + * Orbital Object ToolKit is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; + * without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. + * See the GNU Affero General Public License for more details. + * + * You should have received a copy of the GNU Affero General Public License along with + * Orbital Object ToolKit. If not, see . + */ + +export * from './GroundTrack.js'; +export * from './AccessStatistics.js'; +export * from './AreaCoverage.js'; diff --git a/src/link_budget/LinkBudget.ts b/src/link_budget/LinkBudget.ts new file mode 100644 index 00000000..ce5fe624 --- /dev/null +++ b/src/link_budget/LinkBudget.ts @@ -0,0 +1,314 @@ +/** + * @author Theodore Kruczek + * @description Orbital Object ToolKit (ootk) is a collection of tools for working + * with satellites and other orbital objects. + * @license AGPL-3.0-or-later + * @copyright (c) 2025 Kruczek Labs LLC + * + * Orbital Object ToolKit is free software: you can redistribute it and/or modify it under the + * terms of the GNU Affero General Public License as published by the Free Software + * Foundation, either version 3 of the License, or (at your option) any later version. + * + * Orbital Object ToolKit is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; + * without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. + * See the GNU Affero General Public License for more details. + * + * You should have received a copy of the GNU Affero General Public License along with + * Orbital Object ToolKit. If not, see . + */ + +import { Degrees, Kilometers } from '../types/types.js'; + +/** + * Units for link budget calculations (all in dB unless otherwise specified). + */ +export type Decibels = number; +export type DecibelWatts = number; +export type DecibelMilliwatts = number; +export type DecibelHertz = number; +export type Hertz = number; +export type Watts = number; + +/** + * Transmitter parameters for link budget. + */ +export interface TransmitterParams { + /** Transmit power in dBW */ + powerDbw: DecibelWatts; + /** Antenna gain in dBi */ + antennaGain: Decibels; + /** Line losses in dB (positive value) */ + lineLosses?: Decibels; + /** Pointing loss in dB (positive value) */ + pointingLoss?: Decibels; +} + +/** + * Receiver parameters for link budget. + */ +export interface ReceiverParams { + /** Antenna gain in dBi */ + antennaGain: Decibels; + /** System noise temperature in Kelvin */ + systemNoiseTemp: number; + /** Line losses in dB (positive value) */ + lineLosses?: Decibels; + /** Pointing loss in dB (positive value) */ + pointingLoss?: Decibels; +} + +/** + * Link parameters. + */ +export interface LinkParams { + /** Frequency in Hz */ + frequency: Hertz; + /** Range/distance in kilometers */ + range: Kilometers; + /** Data rate in bits per second */ + dataRate?: number; + /** Atmospheric attenuation in dB (positive value) */ + atmosphericLoss?: Decibels; + /** Rain attenuation in dB (positive value) */ + rainLoss?: Decibels; + /** Ionospheric scintillation loss in dB (positive value) */ + scintillationLoss?: Decibels; + /** Polarization loss in dB (positive value) */ + polarizationLoss?: Decibels; +} + +/** + * Link budget calculation result. + */ +export interface LinkBudgetResult { + /** Effective Isotropic Radiated Power in dBW */ + eirp: DecibelWatts; + /** Free space path loss in dB */ + pathLoss: Decibels; + /** Received power in dBW */ + receivedPower: DecibelWatts; + /** System noise power in dBW */ + noisePower: DecibelWatts; + /** Carrier-to-Noise ratio in dB */ + cnr: Decibels; + /** Carrier-to-Noise density ratio in dB-Hz */ + cn0: DecibelHertz; + /** Energy per bit to noise density ratio in dB (if data rate provided) */ + ebN0?: Decibels; + /** Signal-to-Noise ratio in dB */ + snr: Decibels; + /** Link margin in dB (above minimum required SNR) */ + linkMargin?: Decibels; + /** Total losses in dB */ + totalLosses: Decibels; +} + +/** + * Link budget calculator for RF communication links. + */ +export class LinkBudget { + private static readonly BOLTZMANN_CONSTANT = 1.380649e-23; // J/K + private static readonly SPEED_OF_LIGHT = 299792458; // m/s + + /** + * Calculate link budget for a communication link. + * @param transmitter - Transmitter parameters + * @param receiver - Receiver parameters + * @param link - Link parameters + * @param requiredSnr - Required SNR in dB (optional, for margin calculation) + * @returns Link budget calculation result + */ + static calculate( + transmitter: TransmitterParams, + receiver: ReceiverParams, + link: LinkParams, + requiredSnr?: Decibels, + ): LinkBudgetResult { + // Calculate EIRP + const txLineLoss = transmitter.lineLosses ?? 0; + const txPointingLoss = transmitter.pointingLoss ?? 0; + const eirp = transmitter.powerDbw + transmitter.antennaGain - txLineLoss - txPointingLoss; + + // Calculate free space path loss + const pathLoss = this.calculatePathLoss(link.frequency, link.range); + + // Calculate atmospheric and other losses + const atmosphericLoss = link.atmosphericLoss ?? 0; + const rainLoss = link.rainLoss ?? 0; + const scintillationLoss = link.scintillationLoss ?? 0; + const polarizationLoss = link.polarizationLoss ?? 0; + const rxLineLoss = receiver.lineLosses ?? 0; + const rxPointingLoss = receiver.pointingLoss ?? 0; + + const totalLosses = + pathLoss + atmosphericLoss + rainLoss + scintillationLoss + polarizationLoss + rxLineLoss + rxPointingLoss; + + // Calculate received power + const receivedPower = eirp + receiver.antennaGain - totalLosses; + + // Calculate noise power + const noisePower = this.calculateNoisePower(receiver.systemNoiseTemp, link.dataRate); + + // Calculate CNR + const cnr = receivedPower - noisePower; + + // Calculate C/N0 (carrier to noise density) + const bandwidth = link.dataRate ?? 1; // If no data rate, assume 1 Hz + const cn0 = cnr + 10 * Math.log10(bandwidth); + + // Calculate SNR (assuming matched filter) + const snr = cnr; + + // Calculate Eb/N0 if data rate is provided + let ebN0: Decibels | undefined; + + if (link.dataRate) { + ebN0 = cn0 - 10 * Math.log10(link.dataRate); + } + + // Calculate link margin if required SNR is provided + let linkMargin: Decibels | undefined; + + if (requiredSnr !== undefined) { + linkMargin = snr - requiredSnr; + } + + return { + eirp, + pathLoss, + receivedPower, + noisePower, + cnr, + cn0, + ebN0, + snr, + linkMargin, + totalLosses, + }; + } + + /** + * Calculate free space path loss. + * @param frequency - Frequency in Hz + * @param range - Range in kilometers + * @returns Path loss in dB + */ + static calculatePathLoss(frequency: Hertz, range: Kilometers): Decibels { + const wavelength = this.SPEED_OF_LIGHT / frequency; + const rangeMeters = range * 1000; + + // Friis transmission equation: FSPL = (4πd/λ)² + const pathLoss = 20 * Math.log10((4 * Math.PI * rangeMeters) / wavelength); + + return pathLoss; + } + + /** + * Calculate noise power. + * @param systemNoiseTemp - System noise temperature in Kelvin + * @param bandwidth - Bandwidth in Hz (optional) + * @returns Noise power in dBW + */ + static calculateNoisePower(systemNoiseTemp: number, bandwidth?: number): DecibelWatts { + const bw = bandwidth ?? 1; // If no bandwidth, return noise density + const noisePowerWatts = this.BOLTZMANN_CONSTANT * systemNoiseTemp * bw; + const noisePowerDbw = 10 * Math.log10(noisePowerWatts); + + return noisePowerDbw; + } + + /** + * Calculate antenna gain for a parabolic dish. + * @param diameter - Antenna diameter in meters + * @param frequency - Frequency in Hz + * @param efficiency - Antenna efficiency (0-1, default 0.6) + * @returns Antenna gain in dBi + */ + static calculateParabolicGain(diameter: number, frequency: Hertz, efficiency: number = 0.6): Decibels { + const wavelength = this.SPEED_OF_LIGHT / frequency; + const gain = efficiency * ((Math.PI * diameter) / wavelength) ** 2; + const gainDbi = 10 * Math.log10(gain); + + return gainDbi; + } + + /** + * Calculate G/T (figure of merit for receiver). + * @param antennaGain - Antenna gain in dBi + * @param systemNoiseTemp - System noise temperature in Kelvin + * @returns G/T in dB/K + */ + static calculateGOverT(antennaGain: Decibels, systemNoiseTemp: number): Decibels { + return antennaGain - 10 * Math.log10(systemNoiseTemp); + } + + /** + * Calculate atmospheric attenuation. + * @param frequency - Frequency in Hz + * @param elevationAngle - Elevation angle in degrees + * @param humidity - Relative humidity (0-1, default 0.5) + * @returns Atmospheric attenuation in dB + */ + static calculateAtmosphericLoss(frequency: Hertz, elevationAngle: Degrees, humidity: number = 0.5): Decibels { + // Simplified ITU-R P.676 model for clear air attenuation + const frequencyGHz = frequency / 1e9; + + // Specific attenuation (dB/km) - simplified model + let gamma = 0; + + if (frequencyGHz < 10) { + gamma = 0.01; // Very low below 10 GHz + } else if (frequencyGHz < 20) { + gamma = 0.02 + 0.005 * (frequencyGHz - 10); + } else if (frequencyGHz < 30) { + gamma = 0.07 + 0.01 * (frequencyGHz - 20); + } else { + gamma = 0.17 + 0.02 * Math.min(frequencyGHz - 30, 20); + } + + // Add water vapor contribution + gamma += 0.05 * humidity * (frequencyGHz / 10); + + // Path length through atmosphere (simplified) + const elRad = (elevationAngle * Math.PI) / 180; + const pathLength = 10 / Math.sin(elRad); // Assume 10 km atmosphere height + + return gamma * pathLength; + } + + /** + * Convert dBW to Watts. + * @param dbw - Power in dBW + * @returns Power in Watts + */ + static dbwToWatts(dbw: DecibelWatts): Watts { + return 10 ** (dbw / 10); + } + + /** + * Convert Watts to dBW. + * @param watts - Power in Watts + * @returns Power in dBW + */ + static wattsToDbw(watts: Watts): DecibelWatts { + return 10 * Math.log10(watts); + } + + /** + * Convert dBm to dBW. + * @param dbm - Power in dBm + * @returns Power in dBW + */ + static dbmToDbw(dbm: DecibelMilliwatts): DecibelWatts { + return dbm - 30; + } + + /** + * Convert dBW to dBm. + * @param dbw - Power in dBW + * @returns Power in dBm + */ + static dbwToDbm(dbw: DecibelWatts): DecibelMilliwatts { + return dbw + 30; + } +} diff --git a/src/link_budget/index.ts b/src/link_budget/index.ts new file mode 100644 index 00000000..2f8c5e1c --- /dev/null +++ b/src/link_budget/index.ts @@ -0,0 +1,20 @@ +/** + * @author Theodore Kruczek + * @description Orbital Object ToolKit (ootk) is a collection of tools for working + * with satellites and other orbital objects. + * @license AGPL-3.0-or-later + * @copyright (c) 2025 Kruczek Labs LLC + * + * Orbital Object ToolKit is free software: you can redistribute it and/or modify it under the + * terms of the GNU Affero General Public License as published by the Free Software + * Foundation, either version 3 of the License, or (at your option) any later version. + * + * Orbital Object ToolKit is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; + * without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. + * See the GNU Affero General Public License for more details. + * + * You should have received a copy of the GNU Affero General Public License along with + * Orbital Object ToolKit. If not, see . + */ + +export * from './LinkBudget.js'; diff --git a/src/main.ts b/src/main.ts index a0a8fb48..aa3e7277 100644 --- a/src/main.ts +++ b/src/main.ts @@ -60,3 +60,13 @@ export * from './propagator/index.js'; export * from './orbit_determination/index.js'; export * from './covariance/index.js'; + +export * from './coverage/index.js'; + +export * from './conjunction/index.js'; + +export * from './link_budget/index.js'; + +export * from './constellation/index.js'; + +export * from './attitude/index.js';