Building Blocks of a Timbre Simulation

The example below represents a situation where one has a 21500 second hot dwell at 90.2 degrees pitch that they want to include in a schedule near January 1st, 2021. They want to know how much time at 148.95 degrees pitch that should be included in the schedule to balance this hot dwell from an ACA perspective.

Package Imports

from cxotime import CxoTime
from timbre import get_local_model, find_second_dwell

Define Model Information

This information defines the model being used, as well as initial starting conditions. These starting conditions do not significantly impact the results, however it is good practice to set them to be close to realistic values.

msid = "aacccdpt"
limit = -6.5
date = "2021:001:00:00:00"
aca_model_spec, aca_md5 = get_local_model("../timbre/tests/data/aca_spec.json")
model_init = {"aacccdpt": -6.5, "aca0": -6.5, "eclipse": False}

Define Dwell #1 Information

Define the state information for the observation that one wants to balance, such as a known hot dwell to be included in a schedule. For the ACA model, we only need to define dwell time and pitch. Other models may need to have more information defined. For example, the DPA model would need FEP count, CCD count, SIM position, clocking state, and vid board state defined. Any information that the Xija model needs to run, that isn’t defined in model_init needs to be defined in each of the dwell state definitions, dwell1_state and dwell2_state.

t_dwell1 = 21500.0  # Seconds
dwell1_state = {"pitch": 90.2}

Define Dwell #2 Information

This state can be considered a candidate balancing dwell, such as a cooling attitude one wants to use to balance a hot attitude. In this case we want to calculate the cold time necessary to balance the hot dwell 1 state, so we do not define time as in input parameter. As with dwell1_state, dwell2_state still needs all information, not already in the initialization object model_init, necessary to run the Xija model.

dwell2_state = {"pitch": 148.95}

Calculate Dwell #2 Time

results = find_second_dwell(
    date,
    dwell1_state,
    dwell2_state,
    t_dwell1,
    msid,
    limit,
    aca_model_spec,
    model_init,
    limit_type="max",
)

Explanation of Results

The following information is returned by find_second_dwell:

• converged: This is a boolean value that indicates whether or not a solution was possible. Solutions will not be possible in a number of situations:

  • Both dwells heat the location being modeled.

  • Both dwells cool the location being modeled.

  • One of the states neither sufficiently heats or cools the location being modeled (this will sometimes converge but not always reliably).

  • The dwell1_state is hot and the initial dwell time is long enough to heat this location from a steady state cold temperature to the specified hot limit, assuming the location is associated with a maximum temperature limit. In the case where a location is associated with a minimum temperature limit (e.g. PLINE03T), the opposite would apply, the dwell1_state is cold and is long enough to cool this location from a steady state hot temperature to the specified cold limit.

• unconverged_hot: If the solution didn’t converge, this will be True if the input values resulted in temperatures outside (e.g. above) the specified limit.

• unconverged_cold: If the solution didn’t converge, this will be True if the input values resulted in all temperatures within (e.g. below) the specified limit.

• min_temp: This is the minimum temperature observed during the simulation (latter 2/3 actually to allow the model to reach a repeatable pattern - more on this later). This will be the limit for a converged solution in the case where the location is associated with a minimum temperature limit (e.g. PLINE03T).

• mean_temp: This is the mean temperature observed during the simulation (latter 2/3)

• max_temp: This is the maximum temperature observed during the simulation (latter 2/3). This will be the limit for a converged solution with a maximum temperature limit (e.g. AACCCDPT).

• temperature_limit: This is the limit being used.

• ``dwell2_time``: The dwell #2 time is what you are looking to calculate, and represents the dwell #2 state duration that balances the dwell #1 state at the specified duration, ``t_dwell1``.

• min_pseudo: This is the min pseudo node temperature observed during the evaluated portion of the simulation. This is not implemented yet but is intended to eventually yield additional insight into the results.

• mean_pseudo: This is the mean pseudo node temperature observed during the evaluated portion of the simulation. This is not implemented yet.

• max_pseudo: This is the max pseudo node temperature observed during the evaluated portion of the simulation. This is not implemented yet.

• hotter_state: This is an integer indicating which state is hotter, 1 or 2.

• colder_state: This is an integer indicating which state is colder, 1 or 2, and is actually redundant with hotter_state.

results

{'converged': True,
 'unconverged_hot': False,
 'unconverged_cold': False,
 'min_temp': -7.156585351574264,
 'mean_temp': -6.791883371289597,
 'max_temp': -6.5,
 'temperature_limit': -6.5,
 'dwell_2_time': 47284.94444514031,
 'min_pseudo': nan,
 'mean_pseudo': nan,
 'max_pseudo': nan,
 'hotter_state': 1,
 'colder_state': 2}

Batch Processing

Multiple sets of cases can be run using a single function call, simplifying the generation of larger datasets.

Package Imports

The only difference from above is the import of run_state_pairs and numpy.

import numpy as np
from cxotime import CxoTime
from timbre import get_local_model, run_state_pairs

Define Model Information

This is the same setup as used above.

msid = "aacccdpt"
limit = -6.5
date = "2021:001:00:00:00"
aca_model_spec, aca_md5 = get_local_model("../timbre/tests/data/aca_spec.json")
model_init = {"aacccdpt": -6.5, "aca0": -6.5, "eclipse": False}
t_dwell1 = 21500.0  # Seconds

Define Dwell Cases

This is the most significant departure from above, instead of defining separate dwell1_state and dwell2_state dictionary objects for a single case, pairs of dwell1_state and dwell2_state dictionary objects are combined into a larger data structure. The run_state_pairs function will run through this data structure one pair at a time.

It should be noted that this set of cases will all use the same initial dwell time t_dwell1 listed above.

If one wishes the the return datatype for pitch to be a float, floats should be passed as inputs for all parameters of a given parameter (e.g. pitch).

state_pairs = (
    ({"pitch": 144.2}, {"pitch": 154.95}),
    ({"pitch": 90.2}, {"pitch": 148.95}),
    ({"pitch": 50.0}, {"pitch": 140.0}),
    ({"pitch": 90.0}, {"pitch": 100.0}),
    ({"pitch": 75.0}, {"pitch": 130.0}),
    ({"pitch": 170.0}, {"pitch": 90.0}),
    ({"pitch": 90.0}, {"pitch": 170.0}),
)

Calculate Results

results = run_state_pairs(
    msid,
    aca_model_spec,
    model_init,
    limit,
    date,
    t_dwell1,
    state_pairs,
    limit_type="max",
)

Running simulations for state pair #: 1 out of 7.0

Explanation of Results Format

The description of the results shown above is still valid for the similarly named items, however that information is still included here for completeness, along with descriptions of the additional included information.

• msid: This is the MSID that represents the location of interest, and is the primary output of a given model.

• date: This is the date for which the simulation is applicable.

• datesecs: This is the same date described by date only in seconds using the standard Ska Chandra epoch

• limit: This is the temperature limit being used.

• tdwell_1: This is the initial, “known”, time that corresponds to the dwell #1 state. This is fixed for a batch of state_pairs.

• ``tdwell_2``: The dwell #2 time is what you are looking to calculate, and represents the dwell #2 state duration that balances the dwell #1 state at the specified duration, ``t_dwell1``.

• min_temp: This is the minimum temperature observed during the simulation (latter 2/3 actually to allow the model to reach a repeatable pattern - more on this later). This will be the limit for a converged solution in the case where the location is associated with a minimum temperature limit (e.g. PLINE03T).

• mean_temp: This is the mean temperature observed during the simulation (latter 2/3)

• max_temp: This is the maximum temperature observed during the simulation (latter 2/3). This will be the limit for a converged solution with a maximum temperature limit (e.g. AACCCDPT).

• min_pseudo: This is the min pseudo node temperature observed during the evaluated portion of the simulation. This is not implemented yet but is intended to eventually yield additional insight into the results.

• mean_pseudo: This is the mean pseudo node temperature observed during the evaluated portion of the simulation. This is not implemented yet.

• max_pseudo: This is the max pseudo node temperature observed during the evaluated portion of the simulation. This is not implemented yet.

• converged: This is a boolean value that indicates whether or not a solution was possible. Solutions will not be possible in a number of situations:

  • Both dwells heat the location being modeled.

  • Both dwells cool the location being modeled.

  • One of the states neither sufficiently heats or cools the location being modeled (this will sometimes converge but not always reliably).

  • The dwell1_state is hot and the initial dwell time is long enough to heat this location from a steady state cold temperature to the specified hot limit, assuming the location is associated with a maximum temperature limit. In the case where a location is associated with a minimum temperature limit (e.g. PLINE03T), the opposite would apply, the dwell1_state is cold and is long enough to cool this location from a steady state hot temperature to the specified cold limit.

• unconverged_hot: If the solution didn’t converge, this will be True if the input values resulted in temperatures outside (e.g. above) the specified limit.

• unconverged_cold: If the solution didn’t converge, this will be True if the input values resulted in all temperatures within (e.g. below) the specified limit.

• hotter_state: This is an integer indicating which state is hotter, 1 or 2.

• colder_state: This is an integer indicating which state is colder, 1 or 2, and is actually redundant with hotter_state.

• pitch1: This is the pitch used as an input to a given simulation corresponding to the dwell #1 state.

• eclipse1: This is the eclipse state used as an input to a given simulation corresponding to the dwell #1 state.

• pitch2: This is the pitch used as an input to a given simulation corresponding to the dwell #2 state.

• eclipse2: This is the eclipse state used as an input to a given simulation corresponding to the dwell #2 state.

import astropy

astropy.table.Table(results)

Table length=7

msid	date	datesecs	limit	t_dwell1	t_dwell2	min_temp	mean_temp	max_temp	min_pseudo	mean_pseudo	max_pseudo	converged	unconverged_hot	unconverged_cold	hotter_state	colder_state	pitch1	eclipse1	pitch2	eclipse2
str20	str8	float64	float64	float64	float64	float64	float64	float64	float64	float64	float64	bool	bool	bool	int8	int8	float64	bool	float64	bool
aacccdpt	2021:001	725846469.184	-6.5	21500.0	nan	-8.393848298088113	-8.342559970287693	-8.30616855357068	nan	nan	nan	False	False	True	1	2	144.2	False	154.95	False
aacccdpt	2021:001	725846469.184	-6.5	21500.0	47250.373802043854	-7.154068151351936	-6.790040553689484	-6.5	nan	nan	nan	True	False	False	1	2	90.2	False	148.95	False
aacccdpt	2021:001	725846469.184	-6.5	21500.0	nan	-7.411729330402032	-6.956815970804288	-6.785451637474314	nan	nan	nan	False	False	True	2	1	50.0	False	140.0	False
aacccdpt	2021:001	725846469.184	-6.5	21500.0	nan	-0.012929127833614672	0.03814441333926601	0.10389656662600388	nan	nan	nan	False	True	False	1	2	90.0	False	100.0	False
aacccdpt	2021:001	725846469.184	-6.5	21500.0	nan	-3.2286628278301754	-3.155946443627189	-3.0047592290912712	nan	nan	nan	False	True	False	1	2	75.0	False	130.0	False
aacccdpt	2021:001	725846469.184	-6.5	21500.0	41364.30588583145	-7.6912126491687705	-7.1312558604180865	-6.5	nan	nan	nan	True	False	False	2	1	170.0	False	90.0	False
aacccdpt	2021:001	725846469.184	-6.5	21500.0	10243.34267539528	-6.897765606842965	-6.7052851256062445	-6.5	nan	nan	nan	True	False	False	1	2	90.0	False	170.0	False

table = astropy.table.Table(results)
table.write(
    "run_state_pairs_example_output_table.rst", format="ascii.rst", overwrite=True
)

Discussion of Results

Each number corresponds to a row in the above results:

1. This simulation included a dwell #1 pitch of 144.2 degrees, and a dwell #2 pitch of 154.95 degrees. Although the first dwell state is warmer, this state does not sufficiently heat this model at the given date to reach the specified limit, so these two states together result in an unconverged cold simulation.

2. This simulation included a dwell #1 pitch of 90.2 degrees, and a dwell #2 pitch of 148.95 degrees. This solution converged with a dwell #2 duration of approximately 47250 seconds at 148.95 degrees pitch calculated to sufficiently balance 21500 seconds at 90.2 degrees pitch on 2021:001.

3. This simulation included a dwell #1 pitch of 50.0 degrees, and a dwell #2 pitch of 140.0 degrees. This solution did not converge as both pitch values did not sufficiently heat this model at the given date to reach the specified limit.

4. This simulation included a dwell #1 pitch of 90.0 degrees, and a dwell #2 pitch of 100.0 degrees. As both of these states heat this modeled location, no solution is possible.

5. This simulation included a dwell #1 pitch of 75.0 degrees, and a dwell #2 pitch of 130.0 degrees. As both of these states heat this modeled location, no solution is possible.

6. This simulation included a dwell #1 pitch of 170.0 degrees, and a dwell #2 pitch of 90.0 degrees. This solution converged with a dwell #2 duration of approximately 41364 seconds at 90.0 degrees pitch calculated to sufficiently balance 21500 seconds at 170.0 degrees pitch on 2021:001.

7. This simulation included a dwell #1 pitch of 90.0 degrees, and a dwell #2 pitch of 170.0 degrees. This solution converged with a dwell #2 duration of approximately 10243 seconds at 170.0 degrees pitch calculated to sufficiently balance 21500 seconds at 90.0 degrees pitch on 2021:001.

Background: How Timbre Works

The timbre package enables one to calculate the required contiguous cold time to balance a specified hot observation, or alternatively the maximum contiguous hot time yielded by a given cold observation, for a given Xija thermal model.

Timbre accomplishes this by generating fictitious 30 day schedules that alternate between each of the two observations (and only those two observations), with the duration of the dwell time of interest (t_dwell2) being the only changing variable, for a single Xija thermal model. If a reasonable “dwell time of interest” (t_dwell2) can be found that results in the schedule reaching but not exceeding the specified limit for the thermal model being characterized, then this simulation has converged.

Detailed View of Converged Solution

Below we take a look at such a schedule for the first converged solution from the table above.

• The top plot primarily shows AACCCDPT temperature resulting from this Timbre simulation.

• The bottom plot shows which dwell state is active.

• The thin vertical shaded regions also highlight which dwell state is active.

  • The thin vertical lightly shaded regions are each 21500 seconds, corresponding to the length of t_dwell1.

  • The thin unshaded (white) regions in between the shaded regions are each approximately 53870 seconds, corresponding to the length of t_dwell2, which is the output of a Timbre calculation.

• The dark shaded region on the left highlights the duration used to allow the schedule to reach a “steady state” oscillation, and is not used to evaluate convergence.

• The simulation date is in the center of the evaluation period (2021:001).

• Some small amount of instability in a Xija temperature prediction is not unusual.

• Maximum temperatures in the evaluation region do not need to consistently reach the limit, they only need to be close for the ratio of tdwell_1 to t_dwell2 to be sufficiently accurate for the conditions being used.

• Note that the dwell state transitions do not line up with the peaks and valleys. This is due to the “momentum” built into the model and accurately represents the thermal behavior of location being modeled.

import numpy as np
from cxotime import CxoTime
from timbre import get_local_model, run_state_pairs, create_opt_fun
from timbre import calc_binary_schedule
import plotly.io as pio

from plot_code import *

date = "2021:001:00:00:00"
dwell1_state = {"pitch": 90.2}
dwell2_state = {"pitch": 148.95}
t_dwell1 = 21500.0
msid = "aacccdpt"
limit = -6.5
model_spec, aca_md5 = get_local_model("../timbre/tests/data/aca_spec.json")
init = {"aacccdpt": -6.5, "aca0": -6.5, "eclipse": False}
limit_type = "max"
duration = 2592000
t_backoff = 1725000
n_dwells = 10

datesecs = CxoTime(date).secs

# This ensures three "cycles" of the two dwell states, within the portion of the
# schedule used for evaluation. Subtract 1000 sec for extra padding.
max_dwell = (t_backoff - t_dwell1) / 3 - 1000

t_dwell2 = 47250.37380112706

model_results, state_times, state_keys = calc_binary_schedule(
    datesecs,
    dwell1_state,
    dwell2_state,
    t_dwell1,
    t_dwell2,
    msid,
    model_spec,
    init,
    duration=duration,
    t_backoff=t_backoff,
)


tstart = model_results["aacccdpt"].times[0]
tstop = model_results["aacccdpt"].times[-1]
state_data = {"state_times": state_times, "state_keys": state_keys}
plot_data = format_plot_data(
    model_results["aacccdpt"], -6.5, state_data, dwell1_state, dwell2_state
)
shapes = format_shapes(state_data)
shapes.extend(gen_unused_range(tstart, tstop))
annotations = gen_range_annotations(tstart, tstop, -6.3, -6.3 + 0.01)
annotations.extend(gen_limit_annotation("2021:001", -6.51, -6.5, "Celsius"))
annotations.extend(
    gen_shading_annotation(
        "2021:010:12:00:00",
        -5.7,
        f"Pitch: {dwell1_state['pitch']}",
        f"Pitch: {dwell2_state['pitch']}",
    )
)
plot_dict = generate_converged_solution_plot_dict(
    plot_data, shapes, annotations, tstart, tstop
)
pio.show(plot_dict)

The plot above is intended to provide insight into what a balanced solution looks like, by demonstrating how alternating between two “balanced” dwells results in a relatively steady temperature profile, with the peak temperature reaching but not exceeding the limit. In the interest of minimizing unnecessary data generation this entire profile is not returned by Timbre, instead Timbre returns the duration of the second state, some descriptive information about the solution, and any inputs necessary for characterizing the each state (e.g. pitch, ccd_count, etc.) as described in the Basic and Batch Processing sections.

The plot below demonstrates the conceptual output of Timbre:

plot_dict = generate_example_balance_plot_dict(
    t_dwell1, t_dwell2, dwell1_state, dwell2_state
)
pio.show(plot_dict)

Dwell #2 Calculation Algorithm

Now that a more detailed view of a converged solution has been presented, next we dive into how Timbre arrives at this solution. As stated earlier, the only dependant variable is the duration of the second dwell state, referred to as t_dwell2 in the data returned by Timbre. There are a number of potential approaches to determining this parameter, from brute force to widely used optimization techniques, however the most robust and fastest technique for this particular application has been found to be a three pass brute force + interpolation method: 1. The first pass calculates the values for 1.e-6 and a pre-defined maximum dwell (approximately 500Ks). This serves to eliminate unncessary model runs when no viable solution exists. 2. The second pass calculates the results from N logarithmically (base 10) spaced guesses for t_dwell2, between 1.e-6 and the pre-defined maximum dwell. N is a keyword option, with a default value of 10. A logarithmically spaced set of guesses are used to yield more resolution at the smaller time durations, where the more useful (and likely more accurate) prediction data exists. 3. The third pass calculates the results from N linearly spaced guesses for t_dwell2 between the two closest values calculated from the second pass. 4. The results of the third pass are used to interpolate the value of t_dwell2, using the values for the t_dwell2 guesses and the resulting maximum temperatures observed (in the range of data used for evaluation) to determine the optimum value of t_dwell2 at the limit.

To provide a visual of what the second and third passes look like, please see the following plots:

dwell2_range = np.logspace(1.0e-6, 1, n_dwells, endpoint=True) / n_dwells
dwell2_range = (
    max_dwell * \
    (dwell2_range - dwell2_range[0]) / (dwell2_range[-1] - dwell2_range[0])
)

# t_dwell2 = 47250.37380112706
model_data = {}
for t_guess in dwell2_range:
    model_results, state_times, state_keys = calc_binary_schedule(
        datesecs,
        dwell1_state,
        dwell2_state,
        t_dwell1,
        t_guess,
        msid,
        model_spec,
        init,
        duration=duration,
        t_backoff=t_backoff,
    )
    model_data[t_guess] = dict(
        zip(
            ("model_results", "state_times", "state_keys"),
            (model_results, state_times, state_keys),
        )
    )

tstart = model_data[dwell2_range[0]]["model_results"]["aacccdpt"].times[0]
tstop = model_data[dwell2_range[0]]["model_results"]["aacccdpt"].times[-1]

plot_data = format_step_2_plot_data(model_data, limit, tstart, tstop)
plot_dict = generate_step_2_plot_dict(
    plot_data, tstart, tstop, title="Timbre Second Pass Results", units="Celsius"
)
pio.show(plot_dict)

dwell2_range = np.logspace(1.0e-6, 1, n_dwells, endpoint=True) / n_dwells
dwell2_range = (
    max_dwell * \
    (dwell2_range - dwell2_range[0]) / (dwell2_range[-1] - dwell2_range[0])
)
opt_fun = create_opt_fun(
    datesecs,
    dwell1_state,
    dwell2_state,
    t_dwell1,
    msid,
    model_spec,
    init,
    t_backoff,
    duration,
)
output_pass_2 = np.array(
    [opt_fun(t) for t in dwell2_range],
    dtype=[
        ("duration2", np.float64),
        ("max", np.float64),
        ("mean", np.float64),
        ("min", np.float64),
    ],
)
output_sorted = np.sort(output_pass_2, order="max")
ind = np.searchsorted(output_sorted["max"], limit)

# Generate a linearly spaced set of t_dwell2 guesses from the two closest guesses
# from the second pass
t_bound = (output_sorted["duration2"][ind - 1], output_sorted["duration2"][ind])
dwell2_range = np.linspace(np.min(t_bound), np.max(t_bound), n_dwells,
                           endpoint=True)

model_data = {}
for t_guess in dwell2_range:
    model_results, state_times, state_keys = calc_binary_schedule(
        datesecs,
        dwell1_state,
        dwell2_state,
        t_dwell1,
        t_guess,
        msid,
        model_spec,
        init,
        duration=duration,
        t_backoff=t_backoff,
    )
    model_data[t_guess] = dict(
        zip(
            ("model_results", "state_times", "state_keys", "output"),
            (model_results, state_times, state_keys),
        )
    )

tstart = model_data[dwell2_range[0]]["model_results"]["aacccdpt"].times[0]
tstop = model_data[dwell2_range[0]]["model_results"]["aacccdpt"].times[-1]

plot_data = format_step_2_plot_data(model_data, limit, tstart, tstop)
plot_dict = generate_step_2_plot_dict(
    plot_data, tstart, tstop, title="Timbre Third Pass Results", units="Celsius"
)
pio.show(plot_dict)

This plot is rather busy, making it difficult to estimate the likely optimum second dwell duration by eye. Instead lets look at a plot of just the maximum temperature for each run in the duration used for evaluation (see “Converged Timbre Dwell Simulation” plot above).

output_pass_3 = np.array(
    [opt_fun(t) for t in dwell2_range],
    dtype=[
        ("duration2", np.float64),
        ("max", np.float64),
        ("mean", np.float64),
        ("min", np.float64),
    ],
)

plot_dict = generate_step_3_max_temp_plot_dict(
    output_pass_3,
    "Maximum Temperature for Dwell #2 Guesses (Third Pass)",
    t_dwell2,
    units="Celsius",
)
pio.show(plot_dict)

This plot shows that when starting with a fixed initial dwell of 21500 seconds, at a pitch of 90.2 degrees, the ACA model needs 47250 seconds at 148.95 degrees pitch to maintain a balanced profile. See the plot “Converged Timbre Dwell Simulation” plot above for a more detailed view of the schedule that uses this interpolated dwell #2 time.