Dave Jarvis' Repositories

/* Copyright 2024 White Magic Software, Ltd. -- All rights reserved.
 *
 * SPDX-License-Identifier: MIT
 */
class Layer {
  static TEMPERATURES = {
    TROPOSPHERE: [{
      range: [0, 11],
      temperature: (altitude) => 288.15 - 6.5 * altitude
    }],
    STRATOSPHERE: [{
      range: [11, 20],
      temperature: (_) => 216.65
    }, {
      range: [20, 32],
      temperature: (altitude) => 196.65 + altitude
    }, {
      range: [32, 47],
      temperature: (altitude) => 139.05 + 2.8 * altitude
    }],
    MESOSPHERE: [{
      range: [47, 51],
      temperature: (_) => 270.65
    }, {
      range: [51, 71],
      temperature: (altitude) => 413.45 - 2.8 * altitude
    }, {
      range: [71, 84.852],
      temperature: (altitude) => 356.65 - 2.0 * altitude
    }],
    THERMOSPHERE: [{
      range: [84.852, 91],
      temperature: (_) => 186.8673
    }, {
      range: [91, 110],
      temperature: (altitude) => 263.1905 - 76.3232 * Math.sqrt(1 - Math.pow((altitude - 91) / -19.9429, 2))
    }, {
      range: [110, 120],
      temperature: (altitude) => 240 + 12 * (altitude - 110)
    }, {
      range: [120, 1000],
      temperature: (altitude) => {
        const E = (altitude - 120) * (6356.766 + 120) / (6356.766 + altitude);
        return 1000 - 640 * Math.exp(-0.01875 * E);
      }
    }, {
      range: [1000, Infinity],
      temperature: (_) => 1000
    }]
  };

  constructor(name, range, pressure, density) {
    this.name = name;
    this.range = range;
    this.pressure = pressure;
    this.density = density;
  }

  within(altitude) {
    return altitude >= this.range[0] && altitude <= this.range[1];
  }

  temperature(altitude) {
    const temperatures = Layer.TEMPERATURES[this.name];
    let result;

    for (let t of temperatures) {
      if (altitude >= t.range[0] && altitude <= t.range[1]) {
        result = t.temperature(altitude);
        break;
      }
    }

    console.assert(result !== null && !isNaN(result));
    return result;
  }
}

/**
 * Responsible for computing temperature, air pressure, and air density at a
 * given altitude. Based on Robert A. Braeunig's tables and description:
 *
 * http://www.braeunig.us/space/atmmodel.htm
 */
class Atmosphere {
  static EQUATORIAL_RADIUS = 6378137.0; // m
  static FLATTENING = 1.0 / 298.25722356;
  static SPECIFIC_GAS_CONSTANT = 287.053; // J/kg-K

  static LAYERS = [
    new Layer(
      "TROPOSPHERE",
      [0, 11],
      (altitude) => 101325.0 * Math.pow((288.15 / (288.15 - 6.5 * altitude)), (34.1632 / -6.5)),
      (pressure, gas, temperature) => pressure / (gas * temperature)
    ),
    new Layer(
      "STRATOSPHERE",
      [11, 20],
      (altitude) => 22632.06 * Math.exp(-34.1632 * (altitude - 11) / 216.65),
      (pressure, gas, temperature) => pressure / (gas * temperature)
    ),
    new Layer(
      "STRATOSPHERE",
      [20, 32],
      (altitude) => 5474.889 * Math.pow(216.65 / (216.65 + (altitude - 20)), 34.1632),
      (pressure, gas, temperature) => pressure / (gas * temperature)
    ),
    new Layer(
      "STRATOSPHERE",
      [32, 47],
      (altitude) => 868.0187 * Math.pow(228.65 / (228.65 + 2.8 * (altitude - 32)), 34.1632 / 2.8),
      (pressure, gas, temperature) => pressure / (gas * temperature)
    ),
    new Layer(
      "MESOSPHERE",
      [47, 51],
      (altitude) => 110.9063 * Math.exp(-34.1632 * (altitude - 47) / 270.65),
      (pressure, gas, temperature) => pressure / (gas * temperature)
    ),
    new Layer(
      "MESOSPHERE",
      [51, 71],
      (altitude) => 66.93887 * Math.pow(270.65 / (270.65 - 2.8 * (altitude - 51)), 34.1632 / -2.8),
      (pressure, gas, temperature) => pressure / (gas * temperature)
    ),
    new Layer(
      "MESOSPHERE",
      [71, 84.852],
      (altitude) => 3.956420 * Math.pow(214.65 / (214.65 - 2 * (altitude - 71)), 34.1632 / -2),
      (pressure, gas, temperature) => pressure / (gas * temperature)
    ),
    new Layer(
      "THERMOSPHERE",
      [84.852, Infinity],
      (altitude) =>
        Math.exp(-0.0000000422012 *
          Math.pow(altitude, 5) + 0.0000213489 *
          Math.pow(altitude, 4) - 0.00426388 *
          Math.pow(altitude, 3) + 0.421404 *
          Math.pow(altitude, 2) - 20.8270 * altitude + 416.225),
      (altitude) => Math.exp(
        0.000000075691 *
        Math.pow(altitude, 5) - 0.0000376113 *
        Math.pow(altitude, 4) + 0.0074765 *
        Math.pow(altitude, 3) - 0.743012 *
        Math.pow(altitude, 2) + 36.7280 * altitude - 729.346)
    )
  ];

  constructor(latitude) {
    console.assert(typeof latitude === 'number');
    console.assert(latitude !== undefined);

    this.latitude = latitude;
  }

  radius(latitude) {
    const phi = latitude * (Math.PI / 180);
    return Atmosphere.EQUATORIAL_RADIUS * (1 - Atmosphere.FLATTENING * Math.sin(phi) ** 2);
  }

  geopotential(altitude) {
    const radius = this.radius(this.latitude);
    return (altitude * radius) / (radius + altitude);
  }

  layerValue(altitude, callback) {
    console.assert( altitude >= 0 );

    let result;

    const gpAltitude = this.geopotential(altitude / 1000);

    for (let layer of Atmosphere.LAYERS) {
      if (layer.within(gpAltitude)) {
        result = callback(layer, gpAltitude);
        break;
      }
    }

    this.invariant(result);

    return result;
  }

  temperature(altitude) {
    return this.layerValue(altitude, (layer, gpAltitude) => layer.temperature(gpAltitude));
  }

  airPressure(altitude) {
    return this.layerValue(altitude, (layer, gpAltitude) => layer.pressure(gpAltitude));
  }

  airDensity(altitude) {
    return this.layerValue(altitude, (layer, gpAltitude) => {
      const temperature = layer.temperature(gpAltitude);
      const pressure = layer.pressure(gpAltitude);

      return layer.density(pressure, Atmosphere.SPECIFIC_GAS_CONSTANT, temperature);
    });
  }

  invariant(result) {
    console.assert(result !== null && !isNaN(result));
  }
}

export default Atmosphere;