Skip to content

Stokes I Simultaneous Image and Instrument Modeling

In this tutorial, we will create a preliminary reconstruction of the 2017 M87 data on April 6 by simultaneously creating an image and model for the instrument. By instrument model, we mean something akin to self-calibration in traditional VLBI imaging terminology. However, unlike traditional self-cal, we will solve for the gains each time we update the image self-consistently. This allows us to model the correlations between gains and the image.

To get started we load Comrade.

julia
using Comrade


using Pyehtim
using LinearAlgebra
    CondaPkg Found dependencies: /home/runner/.julia/packages/DimensionalData/hv9KC/CondaPkg.toml
    CondaPkg Found dependencies: /home/runner/.julia/packages/CondaPkg/0UqYV/CondaPkg.toml
    CondaPkg Found dependencies: /home/runner/.julia/packages/PythonCall/83z4q/CondaPkg.toml
    CondaPkg Found dependencies: /home/runner/.julia/packages/Pyehtim/bQtHC/CondaPkg.toml
    CondaPkg Resolving changes
             + ehtim (pip)
             + libstdcxx
             + libstdcxx-ng
             + numpy
             + numpy (pip)
             + openssl
             + pandas
             + python
             + setuptools (pip)
             + uv
             + xarray
    CondaPkg Initialising pixi
             │ /home/runner/.julia/artifacts/cefba4912c2b400756d043a2563ef77a0088866b/bin/pixi
             │ init
             │ --format pixi
             └ /home/runner/work/Comrade.jl/Comrade.jl/examples/intermediate/StokesIImaging/.CondaPkg
✔ Created /home/runner/work/Comrade.jl/Comrade.jl/examples/intermediate/StokesIImaging/.CondaPkg/pixi.toml
    CondaPkg Wrote /home/runner/work/Comrade.jl/Comrade.jl/examples/intermediate/StokesIImaging/.CondaPkg/pixi.toml
             │ [dependencies]
             │ openssl = ">=3, <3.6, >=3, <3.6"
             │ libstdcxx = ">=3.4,<15.0"
             │ uv = ">=0.4"
             │ libstdcxx-ng = ">=3.4,<15.0"
             │ pandas = "<2"
             │ xarray = "*"
             │ numpy = ">=1.24, <2.0"

             │     [dependencies.python]
             │     channel = "conda-forge"
             │     build = "*cp*"
             │     version = ">=3.10,!=3.14.0,!=3.14.1,<4, >=3.6,<=3.12"

             │ [project]
             │ name = ".CondaPkg"
             │ platforms = ["linux-64"]
             │ channels = ["conda-forge"]
             │ channel-priority = "strict"
             │ description = "automatically generated by CondaPkg.jl"

             │ [pypi-dependencies]
             │ ehtim = ">=1.2.10, <2.0"
             │ numpy = ">=1.24, <2.0"
             └ setuptools = "*"
    CondaPkg Installing packages
             │ /home/runner/.julia/artifacts/cefba4912c2b400756d043a2563ef77a0088866b/bin/pixi
             │ install
             └ --manifest-path /home/runner/work/Comrade.jl/Comrade.jl/examples/intermediate/StokesIImaging/.CondaPkg/pixi.toml
✔ The default environment has been installed.
/home/runner/work/Comrade.jl/Comrade.jl/examples/intermediate/StokesIImaging/.CondaPkg/.pixi/envs/default/lib/python3.11/site-packages/ehtim/__init__.py:58: UserWarning: pkg_resources is deprecated as an API. See https://setuptools.pypa.io/en/latest/pkg_resources.html. The pkg_resources package is slated for removal as early as 2025-11-30. Refrain from using this package or pin to Setuptools<81.
  import pkg_resources

For reproducibility we use a stable random number genreator

julia
using StableRNGs
rng = StableRNG(42)
StableRNGs.LehmerRNG(state=0x00000000000000000000000000000055)

Load the Data

To download the data visit https://doi.org/10.25739/g85n-f134 First we will load our data:

julia
obs = ehtim.obsdata.load_uvfits(joinpath(__DIR, "..", "..", "Data", "SR1_M87_2017_096_lo_hops_netcal_StokesI.uvfits"))
Python: <ehtim.obsdata.Obsdata object at 0x7fdf8dcb8850>

Now we do some minor preprocessing:

  • Scan average the data since the data have been preprocessed so that the gain phases coherent.

  • Add 1% systematic noise to deal with calibration issues that cause 1% non-closing errors.

julia
obs = scan_average(obs).add_fractional_noise(0.02)
Python: <ehtim.obsdata.Obsdata object at 0x7fdf87099a50>

Now we extract our complex visibilities.

julia
dvis = extract_table(obs, Visibilities())
EHTObservationTable{Comrade.EHTVisibilityDatum{:I}}
  source:      M87
  mjd:         57849
  bandwidth:   1.856e9
  sites:       [:AA, :AP, :AZ, :JC, :LM, :PV, :SM]
  nsamples:    274

##Building the Model/Posterior

Now, we must build our intensity/visibility model. That is, the model that takes in a named tuple of parameters and perhaps some metadata required to construct the model. For our model, we will use a raster or ContinuousImage for our image model. The model is given below:

The model construction is very similar to Imaging a Black Hole using only Closure Quantities, except we include a large scale gaussian since we want to model the zero baselines. For more information about the image model please read the closure-only example.

julia
function sky(θ, metadata)
    (; fg, c, σimg) = θ
    (; ftot, mimg) = metadata
    # Apply the GMRF fluctuations to the image
    rast = apply_fluctuations(CenteredLR(), mimg, σimg .* c.params)
    pimg = parent(rast)
    @. pimg = (ftot * (1 - fg)) * pimg
    m = ContinuousImage(rast, BSplinePulse{3}())
    # Add a large-scale gaussian to deal with the over-resolved mas flux
    g = modify(Gaussian(), Stretch(μas2rad(500.0), μas2rad(500.0)), Renormalize(ftot * fg))
    x, y = centroid(m)
    return shifted(m, -x, -y) + g
end
sky (generic function with 1 method)

Now, let's set up our image model. The EHT's nominal resolution is 20-25 μas. Additionally, the EHT is not very sensitive to a larger field of view. Typically 60-80 μas is enough to describe the compact flux of M87. Given this, we only need to use a small number of pixels to describe our image.

julia
npix = 48
fovx = μas2rad(200.0)
fovy = μas2rad(200.0)
9.69627362219072e-10

Now let's form our cache's. First, we have our usual image cache which is needed to numerically compute the visibilities.

julia
grid = imagepixels(fovx, fovy, npix, npix)
RectiGrid(
executor: ComradeBase.Serial()
Dimensions: 
(↓ X Sampled{Float64} LinRange{Float64}(-4.747133960864206e-10, 4.747133960864206e-10, 48) ForwardOrdered Regular Points,
→ Y Sampled{Float64} LinRange{Float64}(-4.747133960864206e-10, 4.747133960864206e-10, 48) ForwardOrdered Regular Points)
)

Now we need to specify our image prior. For this work we will use a Gaussian Markov Random field prior Since we are using a Gaussian Markov random field prior we need to first specify our mean image. This behaves somewhat similary to a entropy regularizer in that it will start with an initial guess for the image structure. For this tutorial we will use a a symmetric Gaussian with a FWHM of 50 μas

julia
using VLBIImagePriors
using Distributions
fwhmfac = 2 * sqrt(2 * log(2))
mpr = modify(Gaussian(), Stretch(μas2rad(60.0) ./ fwhmfac))
mimg = intensitymap(mpr, grid)
┌ 48×48 IntensityMap{Float64, 2} ┐
├────────────────────────────────┴─────────────────────────────────────── dims ┐
  ↓ X Sampled{Float64} LinRange{Float64}(-4.747133960864206e-10, 4.747133960864206e-10, 48) ForwardOrdered Regular Points,
  → Y Sampled{Float64} LinRange{Float64}(-4.747133960864206e-10, 4.747133960864206e-10, 48) ForwardOrdered Regular Points
└──────────────────────────────────────────────────────────────────────────────┘
  ↓ →          -4.74713e-10  -4.54513e-10  -4.34312e-10  -4.14112e-10  -3.93911e-10  -3.73711e-10  -3.5351e-10  -3.33309e-10  -3.13109e-10  -2.92908e-10  -2.72708e-10  -2.52507e-10  -2.32307e-10  -2.12106e-10  -1.91905e-10  -1.71705e-10  -1.51504e-10  -1.31304e-10  -1.11103e-10  -9.09026e-11  -7.0702e-11   -5.05014e-11  -3.03009e-11  -1.01003e-11  1.01003e-11  3.03009e-11  5.05014e-11  7.0702e-11   9.09026e-11  1.11103e-10  1.31304e-10  1.51504e-10  1.71705e-10  1.91905e-10  2.12106e-10  2.32307e-10  2.52507e-10  2.72708e-10  2.92908e-10  3.13109e-10  3.33309e-10  3.5351e-10  3.73711e-10  3.93911e-10  4.14112e-10  4.34312e-10  4.54513e-10  4.74713e-10
 -4.74713e-10   1.64166e-9    3.03669e-9    5.46896e-9    9.58947e-9    1.63708e-8    2.72103e-8    4.40334e-8   6.93772e-8    1.06424e-7    1.58944e-7    2.3112e-7     3.27203e-7    4.51007e-7    6.05251e-7    7.90813e-7    1.006e-6      1.24597e-6    1.50247e-6    1.76396e-6    2.01631e-6    2.24395e-6    2.43138e-6    2.56496e-6    2.63448e-6   2.63448e-6   2.56496e-6   2.43138e-6   2.24395e-6   2.01631e-6   1.76396e-6   1.50247e-6   1.24597e-6   1.006e-6     7.90813e-7   6.05251e-7   4.51007e-7   3.27203e-7   2.3112e-7    1.58944e-7   1.06424e-7   6.93772e-8   4.40334e-8  2.72103e-8   1.63708e-8   9.58947e-9   5.46896e-9   3.03669e-9   1.64166e-9
 -4.54513e-10   3.03669e-9    5.61718e-9    1.01163e-8    1.77383e-8    3.02823e-8    5.03328e-8    8.14516e-8   1.28332e-7    1.96859e-7    2.9401e-7     4.2752e-7     6.05251e-7    8.3426e-7     1.11957e-6    1.46282e-6    1.86087e-6    2.30476e-6    2.77922e-6    3.26291e-6    3.72971e-6    4.15078e-6    4.4975e-6     4.74459e-6    4.87318e-6   4.87318e-6   4.74459e-6   4.4975e-6    4.15078e-6   3.72971e-6   3.26291e-6   2.77922e-6   2.30476e-6   1.86087e-6   1.46282e-6   1.11957e-6   8.3426e-7    6.05251e-7   4.2752e-7    2.9401e-7    1.96859e-7   1.28332e-7   8.14516e-8  5.03328e-8   3.02823e-8   1.77383e-8   1.01163e-8   5.61718e-9   3.03669e-9
 -4.34312e-10   5.46896e-9    1.01163e-8    1.82191e-8    3.1946e-8     5.45372e-8    9.06473e-8    1.46691e-7   2.3112e-7     3.54535e-7    5.29501e-7    7.69945e-7    1.09003e-6    1.50247e-6    2.01631e-6    2.63448e-6    3.35135e-6    4.15078e-6    5.00526e-6    5.87637e-6    6.71705e-6    7.47539e-6    8.09982e-6    8.54482e-6    8.77641e-6   8.77641e-6   8.54482e-6   8.09982e-6   7.47539e-6   6.71705e-6   5.87637e-6   5.00526e-6   4.15078e-6   3.35135e-6   2.63448e-6   2.01631e-6   1.50247e-6   1.09003e-6   7.69945e-7   5.29501e-7   3.54535e-7   2.3112e-7    1.46691e-7  9.06473e-8   5.45372e-8   3.1946e-8    1.82191e-8   1.01163e-8   5.46896e-9
 -4.14112e-10   9.58947e-9    1.77383e-8    3.1946e-8     5.60153e-8    9.56274e-8    1.58944e-7    2.57213e-7   4.05255e-7    6.21655e-7    9.28446e-7    1.35005e-6    1.9113e-6     2.63448e-6    3.53547e-6    4.6194e-6     5.87637e-6    7.27813e-6    8.77641e-6    1.03039e-5    1.17779e-5    1.31076e-5    1.42025e-5    1.49828e-5    1.53889e-5   1.53889e-5   1.49828e-5   1.42025e-5   1.31076e-5   1.17779e-5   1.03039e-5   8.77641e-6   7.27813e-6   5.87637e-6   4.6194e-6    3.53547e-6   2.63448e-6   1.9113e-6    1.35005e-6   9.28446e-7   6.21655e-7   4.05255e-7   2.57213e-7  1.58944e-7   9.56274e-8   5.60153e-8   3.1946e-8    1.77383e-8   9.58947e-9
 -3.93911e-10   1.63708e-8    3.02823e-8    5.45372e-8    9.56274e-8    1.63252e-7    2.71345e-7    4.39106e-7   6.91838e-7    1.06127e-6    1.58501e-6    2.30476e-6    3.26291e-6    4.4975e-6     6.03564e-6    7.88608e-6    1.0032e-5     1.2425e-5     1.49828e-5    1.75904e-5    2.01069e-5    2.23769e-5    2.42461e-5    2.55782e-5    2.62714e-5   2.62714e-5   2.55782e-5   2.42461e-5   2.23769e-5   2.01069e-5   1.75904e-5   1.49828e-5   1.2425e-5    1.0032e-5    7.88608e-6   6.03564e-6   4.4975e-6    3.26291e-6   2.30476e-6   1.58501e-6   1.06127e-6   6.91838e-7   4.39106e-7  2.71345e-7   1.63252e-7   9.56274e-8   5.45372e-8   3.02823e-8   1.63708e-8
 -3.73711e-10   2.72103e-8    5.03328e-8    9.06473e-8    1.58944e-7    2.71345e-7    4.51007e-7    7.29848e-7   1.14992e-6    1.76396e-6    2.63448e-6    3.83079e-6    5.42336e-6    7.47539e-6    1.0032e-5     1.31076e-5    1.66743e-5    2.06518e-5    2.49032e-5    2.92374e-5    3.34201e-5    3.71931e-5    4.02999e-5    4.2514e-5     4.36662e-5   4.36662e-5   4.2514e-5    4.02999e-5   3.71931e-5   3.34201e-5   2.92374e-5   2.49032e-5   2.06518e-5   1.66743e-5   1.31076e-5   1.0032e-5    7.47539e-6   5.42336e-6   3.83079e-6   2.63448e-6   1.76396e-6   1.14992e-6   7.29848e-7  4.51007e-7   2.71345e-7   1.58944e-7   9.06473e-8   5.03328e-8   2.72103e-8
 -3.5351e-10    4.40334e-8    8.14516e-8    1.46691e-7    2.57213e-7    4.39106e-7    7.29848e-7    1.18108e-6   1.86087e-6    2.85454e-6    4.26328e-6    6.19922e-6    8.77641e-6    1.20971e-5    1.62343e-5    2.12116e-5    2.69834e-5    3.34201e-5    4.02999e-5    4.73137e-5    5.40824e-5    6.01882e-5    6.52158e-5    6.87987e-5    7.06633e-5   7.06633e-5   6.87987e-5   6.52158e-5   6.01882e-5   5.40824e-5   4.73137e-5   4.02999e-5   3.34201e-5   2.69834e-5   2.12116e-5   1.62343e-5   1.20971e-5   8.77641e-6   6.19922e-6   4.26328e-6   2.85454e-6   1.86087e-6   1.18108e-6  7.29848e-7   4.39106e-7   2.57213e-7   1.46691e-7   8.14516e-8   4.40334e-8
 -3.33309e-10   6.93772e-8    1.28332e-7    2.3112e-7     4.05255e-7    6.91838e-7    1.14992e-6    1.86087e-6   2.93191e-6    4.4975e-6     6.71705e-6    9.76724e-6    1.38277e-5    1.90598e-5    2.55782e-5    3.34201e-5    4.2514e-5     5.26553e-5    6.34949e-5    7.45455e-5    8.521e-5      9.483e-5      0.000102751   0.000108396   0.000111334  0.000111334  0.000108396  0.000102751  9.483e-5     8.521e-5     7.45455e-5   6.34949e-5   5.26553e-5   4.2514e-5    3.34201e-5   2.55782e-5   1.90598e-5   1.38277e-5   9.76724e-6   6.71705e-6   4.4975e-6    2.93191e-6   1.86087e-6  1.14992e-6   6.91838e-7   4.05255e-7   2.3112e-7    1.28332e-7   6.93772e-8
 -3.13109e-10   1.06424e-7    1.96859e-7    3.54535e-7    6.21655e-7    1.06127e-6    1.76396e-6    2.85454e-6   4.4975e-6     6.8991e-6     1.03039e-5    1.49828e-5    2.12116e-5    2.92374e-5    3.92365e-5    5.12659e-5    6.52158e-5    8.07724e-5    9.74002e-5    0.000114352   0.000130711   0.000145468   0.000157619   0.000166278   0.000170785  0.000170785  0.000166278  0.000157619  0.000145468  0.000130711  0.000114352  9.74002e-5   8.07724e-5   6.52158e-5   5.12659e-5   3.92365e-5   2.92374e-5   2.12116e-5   1.49828e-5   1.03039e-5   6.8991e-6    4.4975e-6    2.85454e-6  1.76396e-6   1.06127e-6   6.21655e-7   3.54535e-7   1.96859e-7   1.06424e-7
 -2.92908e-10   1.58944e-7    2.9401e-7     5.29501e-7    9.28446e-7    1.58501e-6    2.63448e-6    4.26328e-6   6.71705e-6    1.03039e-5    1.53889e-5    2.23769e-5    3.16796e-5    4.36662e-5    5.86e-5       7.65659e-5    9.74002e-5    0.000120634   0.000145468   0.000170785   0.000195218   0.000217257   0.000235405   0.000248338   0.000255069  0.000255069  0.000248338  0.000235405  0.000217257  0.000195218  0.000170785  0.000145468  0.000120634  9.74002e-5   7.65659e-5   5.86e-5      4.36662e-5   3.16796e-5   2.23769e-5   1.53889e-5   1.03039e-5   6.71705e-6   4.26328e-6  2.63448e-6   1.58501e-6   9.28446e-7   5.29501e-7   2.9401e-7    1.58944e-7
 -2.72708e-10   2.3112e-7     4.2752e-7     7.69945e-7    1.35005e-6    2.30476e-6    3.83079e-6    6.19922e-6   9.76724e-6    1.49828e-5    2.23769e-5    3.25382e-5    4.60652e-5    6.34949e-5    8.521e-5      0.000111334   0.000141629   0.000175414   0.000211524   0.000248338   0.000283865   0.000315913   0.000342302   0.000361107   0.000370894  0.000370894  0.000361107  0.000342302  0.000315913  0.000283865  0.000248338  0.000211524  0.000175414  0.000141629  0.000111334  8.521e-5     6.34949e-5   4.60652e-5   3.25382e-5   2.23769e-5   1.49828e-5   9.76724e-6   6.19922e-6  3.83079e-6   2.30476e-6   1.35005e-6   7.69945e-7   4.2752e-7    2.3112e-7
 -2.52507e-10   3.27203e-7    6.05251e-7    1.09003e-6    1.9113e-6     3.26291e-6    5.42336e-6    8.77641e-6   1.38277e-5    2.12116e-5    3.16796e-5    4.60652e-5    6.52158e-5    8.98914e-5    0.000120634   0.000157619   0.000200508   0.000248338   0.000299461   0.000351579   0.000401876   0.000447247   0.000484606   0.00051123    0.000525085  0.000525085  0.00051123   0.000484606  0.000447247  0.000401876  0.000351579  0.000299461  0.000248338  0.000200508  0.000157619  0.000120634  8.98914e-5   6.52158e-5   4.60652e-5   3.16796e-5   2.12116e-5   1.38277e-5   8.77641e-6  5.42336e-6   3.26291e-6   1.9113e-6    1.09003e-6   6.05251e-7   3.27203e-7
 -2.32307e-10   4.51007e-7    8.3426e-7     1.50247e-6    2.63448e-6    4.4975e-6     7.47539e-6    1.20971e-5   1.90598e-5    2.92374e-5    4.36662e-5    6.34949e-5    8.98914e-5    0.000123904   0.000166278   0.000217257   0.000276375   0.000342302   0.000412768   0.000484606   0.000553933   0.000616471   0.000667966   0.000704663   0.000723762  0.000723762  0.000704663  0.000667966  0.000616471  0.000553933  0.000484606  0.000412768  0.000342302  0.000276375  0.000217257  0.000166278  0.000123904  8.98914e-5   6.34949e-5   4.36662e-5   2.92374e-5   1.90598e-5   1.20971e-5  7.47539e-6   4.4975e-6    2.63448e-6   1.50247e-6   8.3426e-7    4.51007e-7
 -2.12106e-10   6.05251e-7    1.11957e-6    2.01631e-6    3.53547e-6    6.03564e-6    1.0032e-5     1.62343e-5   2.55782e-5    3.92365e-5    5.86e-5       8.521e-5      0.000120634   0.000166278   0.000223145   0.000291559   0.000370894   0.000459368   0.000553933   0.00065034    0.000743377   0.000827303   0.000896409   0.000945657   0.000971287  0.000971287  0.000945657  0.000896409  0.000827303  0.000743377  0.00065034   0.000553933  0.000459368  0.000370894  0.000291559  0.000223145  0.000166278  0.000120634  8.521e-5     5.86e-5      3.92365e-5   2.55782e-5   1.62343e-5  1.0032e-5    6.03564e-6   3.53547e-6   2.01631e-6   1.11957e-6   6.05251e-7
 -1.91905e-10   7.90813e-7    1.46282e-6    2.63448e-6    4.6194e-6     7.88608e-6    1.31076e-5    2.12116e-5   3.34201e-5    5.12659e-5    7.65659e-5    0.000111334   0.000157619   0.000217257   0.000291559   0.000380947   0.000484606   0.000600204   0.000723762   0.000849725   0.000971287   0.00108094    0.00117124    0.00123558    0.00126907   0.00126907   0.00123558   0.00117124   0.00108094   0.000971287  0.000849725  0.000723762  0.000600204  0.000484606  0.000380947  0.000291559  0.000217257  0.000157619  0.000111334  7.65659e-5   5.12659e-5   3.34201e-5   2.12116e-5  1.31076e-5   7.88608e-6   4.6194e-6    2.63448e-6   1.46282e-6   7.90813e-7
 -1.71705e-10   1.006e-6      1.86087e-6    3.35135e-6    5.87637e-6    1.0032e-5     1.66743e-5    2.69834e-5   4.2514e-5     6.52158e-5    9.74002e-5    0.000141629   0.000200508   0.000276375   0.000370894   0.000484606   0.000616471   0.000763525   0.000920704   0.00108094    0.00123558    0.00137508    0.00148994    0.0015718     0.0016144    0.0016144    0.0015718    0.00148994   0.00137508   0.00123558   0.00108094   0.000920704  0.000763525  0.000616471  0.000484606  0.000370894  0.000276375  0.000200508  0.000141629  9.74002e-5   6.52158e-5   4.2514e-5    2.69834e-5  1.66743e-5   1.0032e-5    5.87637e-6   3.35135e-6   1.86087e-6   1.006e-6
 -1.51504e-10   1.24597e-6    2.30476e-6    4.15078e-6    7.27813e-6    1.2425e-5     2.06518e-5    3.34201e-5   5.26553e-5    8.07724e-5    0.000120634   0.000175414   0.000248338   0.000342302   0.000459368   0.000600204   0.000763525   0.000945657   0.00114033    0.00133879    0.00153032    0.00170309    0.00184535    0.00194673    0.0019995    0.0019995    0.00194673   0.00184535   0.00170309   0.00153032   0.00133879   0.00114033   0.000945657  0.000763525  0.000600204  0.000459368  0.000342302  0.000248338  0.000175414  0.000120634  8.07724e-5   5.26553e-5   3.34201e-5  2.06518e-5   1.2425e-5    7.27813e-6   4.15078e-6   2.30476e-6   1.24597e-6
 -1.31304e-10   1.50247e-6    2.77922e-6    5.00526e-6    8.77641e-6    1.49828e-5    2.49032e-5    4.02999e-5   6.34949e-5    9.74002e-5    0.000145468   0.000211524   0.000299461   0.000412768   0.000553933   0.000723762   0.000920704   0.00114033    0.00137508    0.0016144     0.00184535    0.00205369    0.00222523    0.00234749    0.00241111   0.00241111   0.00234749   0.00222523   0.00205369   0.00184535   0.0016144    0.00137508   0.00114033   0.000920704  0.000723762  0.000553933  0.000412768  0.000299461  0.000211524  0.000145468  9.74002e-5   6.34949e-5   4.02999e-5  2.49032e-5   1.49828e-5   8.77641e-6   5.00526e-6   2.77922e-6   1.50247e-6
 -1.11103e-10   1.76396e-6    3.26291e-6    5.87637e-6    1.03039e-5    1.75904e-5    2.92374e-5    4.73137e-5   7.45455e-5    0.000114352   0.000170785   0.000248338   0.000351579   0.000484606   0.00065034    0.000849725   0.00108094    0.00133879    0.0016144     0.00189536    0.00216652    0.00241111    0.00261251    0.00275604    0.00283074   0.00283074   0.00275604   0.00261251   0.00241111   0.00216652   0.00189536   0.0016144    0.00133879   0.00108094   0.000849725  0.00065034   0.000484606  0.000351579  0.000248338  0.000170785  0.000114352  7.45455e-5   4.73137e-5  2.92374e-5   1.75904e-5   1.03039e-5   5.87637e-6   3.26291e-6   1.76396e-6
 -9.09026e-11   2.01631e-6    3.72971e-6    6.71705e-6    1.17779e-5    2.01069e-5    3.34201e-5    5.40824e-5   8.521e-5      0.000130711   0.000195218   0.000283865   0.000401876   0.000553933   0.000743377   0.000971287   0.00123558    0.00153032    0.00184535    0.00216652    0.00247646    0.00275604    0.00298626    0.00315032    0.00323571   0.00323571   0.00315032   0.00298626   0.00275604   0.00247646   0.00216652   0.00184535   0.00153032   0.00123558   0.000971287  0.000743377  0.000553933  0.000401876  0.000283865  0.000195218  0.000130711  8.521e-5     5.40824e-5  3.34201e-5   2.01069e-5   1.17779e-5   6.71705e-6   3.72971e-6   2.01631e-6
 -7.0702e-11    2.24395e-6    4.15078e-6    7.47539e-6    1.31076e-5    2.23769e-5    3.71931e-5    6.01882e-5   9.483e-5      0.000145468   0.000217257   0.000315913   0.000447247   0.000616471   0.000827303   0.00108094    0.00137508    0.00170309    0.00205369    0.00241111    0.00275604    0.0030672     0.0033234     0.00350599    0.00360101   0.00360101   0.00350599   0.0033234    0.0030672    0.00275604   0.00241111   0.00205369   0.00170309   0.00137508   0.00108094   0.000827303  0.000616471  0.000447247  0.000315913  0.000217257  0.000145468  9.483e-5     6.01882e-5  3.71931e-5   2.23769e-5   1.31076e-5   7.47539e-6   4.15078e-6   2.24395e-6
 -5.05014e-11   2.43138e-6    4.4975e-6     8.09982e-6    1.42025e-5    2.42461e-5    4.02999e-5    6.52158e-5   0.000102751   0.000157619   0.000235405   0.000342302   0.000484606   0.000667966   0.000896409   0.00117124    0.00148994    0.00184535    0.00222523    0.00261251    0.00298626    0.0033234     0.00360101    0.00379885    0.00390181   0.00390181   0.00379885   0.00360101   0.0033234    0.00298626   0.00261251   0.00222523   0.00184535   0.00148994   0.00117124   0.000896409  0.000667966  0.000484606  0.000342302  0.000235405  0.000157619  0.000102751  6.52158e-5  4.02999e-5   2.42461e-5   1.42025e-5   8.09982e-6   4.4975e-6    2.43138e-6
 -3.03009e-11   2.56496e-6    4.74459e-6    8.54482e-6    1.49828e-5    2.55782e-5    4.2514e-5     6.87987e-5   0.000108396   0.000166278   0.000248338   0.000361107   0.00051123    0.000704663   0.000945657   0.00123558    0.0015718     0.00194673    0.00234749    0.00275604    0.00315032    0.00350599    0.00379885    0.00400756    0.00411617   0.00411617   0.00400756   0.00379885   0.00350599   0.00315032   0.00275604   0.00234749   0.00194673   0.0015718    0.00123558   0.000945657  0.000704663  0.00051123   0.000361107  0.000248338  0.000166278  0.000108396  6.87987e-5  4.2514e-5    2.55782e-5   1.49828e-5   8.54482e-6   4.74459e-6   2.56496e-6
 -1.01003e-11   2.63448e-6    4.87318e-6    8.77641e-6    1.53889e-5    2.62714e-5    4.36662e-5    7.06633e-5   0.000111334   0.000170785   0.000255069   0.000370894   0.000525085   0.000723762   0.000971287   0.00126907    0.0016144     0.0019995     0.00241111    0.00283074    0.00323571    0.00360101    0.00390181    0.00411617    0.00422773   0.00422773   0.00411617   0.00390181   0.00360101   0.00323571   0.00283074   0.00241111   0.0019995    0.0016144    0.00126907   0.000971287  0.000723762  0.000525085  0.000370894  0.000255069  0.000170785  0.000111334  7.06633e-5  4.36662e-5   2.62714e-5   1.53889e-5   8.77641e-6   4.87318e-6   2.63448e-6
  1.01003e-11   2.63448e-6    4.87318e-6    8.77641e-6    1.53889e-5    2.62714e-5    4.36662e-5    7.06633e-5   0.000111334   0.000170785   0.000255069   0.000370894   0.000525085   0.000723762   0.000971287   0.00126907    0.0016144     0.0019995     0.00241111    0.00283074    0.00323571    0.00360101    0.00390181    0.00411617    0.00422773   0.00422773   0.00411617   0.00390181   0.00360101   0.00323571   0.00283074   0.00241111   0.0019995    0.0016144    0.00126907   0.000971287  0.000723762  0.000525085  0.000370894  0.000255069  0.000170785  0.000111334  7.06633e-5  4.36662e-5   2.62714e-5   1.53889e-5   8.77641e-6   4.87318e-6   2.63448e-6
  3.03009e-11   2.56496e-6    4.74459e-6    8.54482e-6    1.49828e-5    2.55782e-5    4.2514e-5     6.87987e-5   0.000108396   0.000166278   0.000248338   0.000361107   0.00051123    0.000704663   0.000945657   0.00123558    0.0015718     0.00194673    0.00234749    0.00275604    0.00315032    0.00350599    0.00379885    0.00400756    0.00411617   0.00411617   0.00400756   0.00379885   0.00350599   0.00315032   0.00275604   0.00234749   0.00194673   0.0015718    0.00123558   0.000945657  0.000704663  0.00051123   0.000361107  0.000248338  0.000166278  0.000108396  6.87987e-5  4.2514e-5    2.55782e-5   1.49828e-5   8.54482e-6   4.74459e-6   2.56496e-6
  5.05014e-11   2.43138e-6    4.4975e-6     8.09982e-6    1.42025e-5    2.42461e-5    4.02999e-5    6.52158e-5   0.000102751   0.000157619   0.000235405   0.000342302   0.000484606   0.000667966   0.000896409   0.00117124    0.00148994    0.00184535    0.00222523    0.00261251    0.00298626    0.0033234     0.00360101    0.00379885    0.00390181   0.00390181   0.00379885   0.00360101   0.0033234    0.00298626   0.00261251   0.00222523   0.00184535   0.00148994   0.00117124   0.000896409  0.000667966  0.000484606  0.000342302  0.000235405  0.000157619  0.000102751  6.52158e-5  4.02999e-5   2.42461e-5   1.42025e-5   8.09982e-6   4.4975e-6    2.43138e-6
  7.0702e-11    2.24395e-6    4.15078e-6    7.47539e-6    1.31076e-5    2.23769e-5    3.71931e-5    6.01882e-5   9.483e-5      0.000145468   0.000217257   0.000315913   0.000447247   0.000616471   0.000827303   0.00108094    0.00137508    0.00170309    0.00205369    0.00241111    0.00275604    0.0030672     0.0033234     0.00350599    0.00360101   0.00360101   0.00350599   0.0033234    0.0030672    0.00275604   0.00241111   0.00205369   0.00170309   0.00137508   0.00108094   0.000827303  0.000616471  0.000447247  0.000315913  0.000217257  0.000145468  9.483e-5     6.01882e-5  3.71931e-5   2.23769e-5   1.31076e-5   7.47539e-6   4.15078e-6   2.24395e-6
  9.09026e-11   2.01631e-6    3.72971e-6    6.71705e-6    1.17779e-5    2.01069e-5    3.34201e-5    5.40824e-5   8.521e-5      0.000130711   0.000195218   0.000283865   0.000401876   0.000553933   0.000743377   0.000971287   0.00123558    0.00153032    0.00184535    0.00216652    0.00247646    0.00275604    0.00298626    0.00315032    0.00323571   0.00323571   0.00315032   0.00298626   0.00275604   0.00247646   0.00216652   0.00184535   0.00153032   0.00123558   0.000971287  0.000743377  0.000553933  0.000401876  0.000283865  0.000195218  0.000130711  8.521e-5     5.40824e-5  3.34201e-5   2.01069e-5   1.17779e-5   6.71705e-6   3.72971e-6   2.01631e-6
  1.11103e-10   1.76396e-6    3.26291e-6    5.87637e-6    1.03039e-5    1.75904e-5    2.92374e-5    4.73137e-5   7.45455e-5    0.000114352   0.000170785   0.000248338   0.000351579   0.000484606   0.00065034    0.000849725   0.00108094    0.00133879    0.0016144     0.00189536    0.00216652    0.00241111    0.00261251    0.00275604    0.00283074   0.00283074   0.00275604   0.00261251   0.00241111   0.00216652   0.00189536   0.0016144    0.00133879   0.00108094   0.000849725  0.00065034   0.000484606  0.000351579  0.000248338  0.000170785  0.000114352  7.45455e-5   4.73137e-5  2.92374e-5   1.75904e-5   1.03039e-5   5.87637e-6   3.26291e-6   1.76396e-6
  1.31304e-10   1.50247e-6    2.77922e-6    5.00526e-6    8.77641e-6    1.49828e-5    2.49032e-5    4.02999e-5   6.34949e-5    9.74002e-5    0.000145468   0.000211524   0.000299461   0.000412768   0.000553933   0.000723762   0.000920704   0.00114033    0.00137508    0.0016144     0.00184535    0.00205369    0.00222523    0.00234749    0.00241111   0.00241111   0.00234749   0.00222523   0.00205369   0.00184535   0.0016144    0.00137508   0.00114033   0.000920704  0.000723762  0.000553933  0.000412768  0.000299461  0.000211524  0.000145468  9.74002e-5   6.34949e-5   4.02999e-5  2.49032e-5   1.49828e-5   8.77641e-6   5.00526e-6   2.77922e-6   1.50247e-6
  1.51504e-10   1.24597e-6    2.30476e-6    4.15078e-6    7.27813e-6    1.2425e-5     2.06518e-5    3.34201e-5   5.26553e-5    8.07724e-5    0.000120634   0.000175414   0.000248338   0.000342302   0.000459368   0.000600204   0.000763525   0.000945657   0.00114033    0.00133879    0.00153032    0.00170309    0.00184535    0.00194673    0.0019995    0.0019995    0.00194673   0.00184535   0.00170309   0.00153032   0.00133879   0.00114033   0.000945657  0.000763525  0.000600204  0.000459368  0.000342302  0.000248338  0.000175414  0.000120634  8.07724e-5   5.26553e-5   3.34201e-5  2.06518e-5   1.2425e-5    7.27813e-6   4.15078e-6   2.30476e-6   1.24597e-6
  1.71705e-10   1.006e-6      1.86087e-6    3.35135e-6    5.87637e-6    1.0032e-5     1.66743e-5    2.69834e-5   4.2514e-5     6.52158e-5    9.74002e-5    0.000141629   0.000200508   0.000276375   0.000370894   0.000484606   0.000616471   0.000763525   0.000920704   0.00108094    0.00123558    0.00137508    0.00148994    0.0015718     0.0016144    0.0016144    0.0015718    0.00148994   0.00137508   0.00123558   0.00108094   0.000920704  0.000763525  0.000616471  0.000484606  0.000370894  0.000276375  0.000200508  0.000141629  9.74002e-5   6.52158e-5   4.2514e-5    2.69834e-5  1.66743e-5   1.0032e-5    5.87637e-6   3.35135e-6   1.86087e-6   1.006e-6
  1.91905e-10   7.90813e-7    1.46282e-6    2.63448e-6    4.6194e-6     7.88608e-6    1.31076e-5    2.12116e-5   3.34201e-5    5.12659e-5    7.65659e-5    0.000111334   0.000157619   0.000217257   0.000291559   0.000380947   0.000484606   0.000600204   0.000723762   0.000849725   0.000971287   0.00108094    0.00117124    0.00123558    0.00126907   0.00126907   0.00123558   0.00117124   0.00108094   0.000971287  0.000849725  0.000723762  0.000600204  0.000484606  0.000380947  0.000291559  0.000217257  0.000157619  0.000111334  7.65659e-5   5.12659e-5   3.34201e-5   2.12116e-5  1.31076e-5   7.88608e-6   4.6194e-6    2.63448e-6   1.46282e-6   7.90813e-7
  2.12106e-10   6.05251e-7    1.11957e-6    2.01631e-6    3.53547e-6    6.03564e-6    1.0032e-5     1.62343e-5   2.55782e-5    3.92365e-5    5.86e-5       8.521e-5      0.000120634   0.000166278   0.000223145   0.000291559   0.000370894   0.000459368   0.000553933   0.00065034    0.000743377   0.000827303   0.000896409   0.000945657   0.000971287  0.000971287  0.000945657  0.000896409  0.000827303  0.000743377  0.00065034   0.000553933  0.000459368  0.000370894  0.000291559  0.000223145  0.000166278  0.000120634  8.521e-5     5.86e-5      3.92365e-5   2.55782e-5   1.62343e-5  1.0032e-5    6.03564e-6   3.53547e-6   2.01631e-6   1.11957e-6   6.05251e-7
  2.32307e-10   4.51007e-7    8.3426e-7     1.50247e-6    2.63448e-6    4.4975e-6     7.47539e-6    1.20971e-5   1.90598e-5    2.92374e-5    4.36662e-5    6.34949e-5    8.98914e-5    0.000123904   0.000166278   0.000217257   0.000276375   0.000342302   0.000412768   0.000484606   0.000553933   0.000616471   0.000667966   0.000704663   0.000723762  0.000723762  0.000704663  0.000667966  0.000616471  0.000553933  0.000484606  0.000412768  0.000342302  0.000276375  0.000217257  0.000166278  0.000123904  8.98914e-5   6.34949e-5   4.36662e-5   2.92374e-5   1.90598e-5   1.20971e-5  7.47539e-6   4.4975e-6    2.63448e-6   1.50247e-6   8.3426e-7    4.51007e-7
  2.52507e-10   3.27203e-7    6.05251e-7    1.09003e-6    1.9113e-6     3.26291e-6    5.42336e-6    8.77641e-6   1.38277e-5    2.12116e-5    3.16796e-5    4.60652e-5    6.52158e-5    8.98914e-5    0.000120634   0.000157619   0.000200508   0.000248338   0.000299461   0.000351579   0.000401876   0.000447247   0.000484606   0.00051123    0.000525085  0.000525085  0.00051123   0.000484606  0.000447247  0.000401876  0.000351579  0.000299461  0.000248338  0.000200508  0.000157619  0.000120634  8.98914e-5   6.52158e-5   4.60652e-5   3.16796e-5   2.12116e-5   1.38277e-5   8.77641e-6  5.42336e-6   3.26291e-6   1.9113e-6    1.09003e-6   6.05251e-7   3.27203e-7
  2.72708e-10   2.3112e-7     4.2752e-7     7.69945e-7    1.35005e-6    2.30476e-6    3.83079e-6    6.19922e-6   9.76724e-6    1.49828e-5    2.23769e-5    3.25382e-5    4.60652e-5    6.34949e-5    8.521e-5      0.000111334   0.000141629   0.000175414   0.000211524   0.000248338   0.000283865   0.000315913   0.000342302   0.000361107   0.000370894  0.000370894  0.000361107  0.000342302  0.000315913  0.000283865  0.000248338  0.000211524  0.000175414  0.000141629  0.000111334  8.521e-5     6.34949e-5   4.60652e-5   3.25382e-5   2.23769e-5   1.49828e-5   9.76724e-6   6.19922e-6  3.83079e-6   2.30476e-6   1.35005e-6   7.69945e-7   4.2752e-7    2.3112e-7
  2.92908e-10   1.58944e-7    2.9401e-7     5.29501e-7    9.28446e-7    1.58501e-6    2.63448e-6    4.26328e-6   6.71705e-6    1.03039e-5    1.53889e-5    2.23769e-5    3.16796e-5    4.36662e-5    5.86e-5       7.65659e-5    9.74002e-5    0.000120634   0.000145468   0.000170785   0.000195218   0.000217257   0.000235405   0.000248338   0.000255069  0.000255069  0.000248338  0.000235405  0.000217257  0.000195218  0.000170785  0.000145468  0.000120634  9.74002e-5   7.65659e-5   5.86e-5      4.36662e-5   3.16796e-5   2.23769e-5   1.53889e-5   1.03039e-5   6.71705e-6   4.26328e-6  2.63448e-6   1.58501e-6   9.28446e-7   5.29501e-7   2.9401e-7    1.58944e-7
  3.13109e-10   1.06424e-7    1.96859e-7    3.54535e-7    6.21655e-7    1.06127e-6    1.76396e-6    2.85454e-6   4.4975e-6     6.8991e-6     1.03039e-5    1.49828e-5    2.12116e-5    2.92374e-5    3.92365e-5    5.12659e-5    6.52158e-5    8.07724e-5    9.74002e-5    0.000114352   0.000130711   0.000145468   0.000157619   0.000166278   0.000170785  0.000170785  0.000166278  0.000157619  0.000145468  0.000130711  0.000114352  9.74002e-5   8.07724e-5   6.52158e-5   5.12659e-5   3.92365e-5   2.92374e-5   2.12116e-5   1.49828e-5   1.03039e-5   6.8991e-6    4.4975e-6    2.85454e-6  1.76396e-6   1.06127e-6   6.21655e-7   3.54535e-7   1.96859e-7   1.06424e-7
  3.33309e-10   6.93772e-8    1.28332e-7    2.3112e-7     4.05255e-7    6.91838e-7    1.14992e-6    1.86087e-6   2.93191e-6    4.4975e-6     6.71705e-6    9.76724e-6    1.38277e-5    1.90598e-5    2.55782e-5    3.34201e-5    4.2514e-5     5.26553e-5    6.34949e-5    7.45455e-5    8.521e-5      9.483e-5      0.000102751   0.000108396   0.000111334  0.000111334  0.000108396  0.000102751  9.483e-5     8.521e-5     7.45455e-5   6.34949e-5   5.26553e-5   4.2514e-5    3.34201e-5   2.55782e-5   1.90598e-5   1.38277e-5   9.76724e-6   6.71705e-6   4.4975e-6    2.93191e-6   1.86087e-6  1.14992e-6   6.91838e-7   4.05255e-7   2.3112e-7    1.28332e-7   6.93772e-8
  3.5351e-10    4.40334e-8    8.14516e-8    1.46691e-7    2.57213e-7    4.39106e-7    7.29848e-7    1.18108e-6   1.86087e-6    2.85454e-6    4.26328e-6    6.19922e-6    8.77641e-6    1.20971e-5    1.62343e-5    2.12116e-5    2.69834e-5    3.34201e-5    4.02999e-5    4.73137e-5    5.40824e-5    6.01882e-5    6.52158e-5    6.87987e-5    7.06633e-5   7.06633e-5   6.87987e-5   6.52158e-5   6.01882e-5   5.40824e-5   4.73137e-5   4.02999e-5   3.34201e-5   2.69834e-5   2.12116e-5   1.62343e-5   1.20971e-5   8.77641e-6   6.19922e-6   4.26328e-6   2.85454e-6   1.86087e-6   1.18108e-6  7.29848e-7   4.39106e-7   2.57213e-7   1.46691e-7   8.14516e-8   4.40334e-8
  3.73711e-10   2.72103e-8    5.03328e-8    9.06473e-8    1.58944e-7    2.71345e-7    4.51007e-7    7.29848e-7   1.14992e-6    1.76396e-6    2.63448e-6    3.83079e-6    5.42336e-6    7.47539e-6    1.0032e-5     1.31076e-5    1.66743e-5    2.06518e-5    2.49032e-5    2.92374e-5    3.34201e-5    3.71931e-5    4.02999e-5    4.2514e-5     4.36662e-5   4.36662e-5   4.2514e-5    4.02999e-5   3.71931e-5   3.34201e-5   2.92374e-5   2.49032e-5   2.06518e-5   1.66743e-5   1.31076e-5   1.0032e-5    7.47539e-6   5.42336e-6   3.83079e-6   2.63448e-6   1.76396e-6   1.14992e-6   7.29848e-7  4.51007e-7   2.71345e-7   1.58944e-7   9.06473e-8   5.03328e-8   2.72103e-8
  3.93911e-10   1.63708e-8    3.02823e-8    5.45372e-8    9.56274e-8    1.63252e-7    2.71345e-7    4.39106e-7   6.91838e-7    1.06127e-6    1.58501e-6    2.30476e-6    3.26291e-6    4.4975e-6     6.03564e-6    7.88608e-6    1.0032e-5     1.2425e-5     1.49828e-5    1.75904e-5    2.01069e-5    2.23769e-5    2.42461e-5    2.55782e-5    2.62714e-5   2.62714e-5   2.55782e-5   2.42461e-5   2.23769e-5   2.01069e-5   1.75904e-5   1.49828e-5   1.2425e-5    1.0032e-5    7.88608e-6   6.03564e-6   4.4975e-6    3.26291e-6   2.30476e-6   1.58501e-6   1.06127e-6   6.91838e-7   4.39106e-7  2.71345e-7   1.63252e-7   9.56274e-8   5.45372e-8   3.02823e-8   1.63708e-8
  4.14112e-10   9.58947e-9    1.77383e-8    3.1946e-8     5.60153e-8    9.56274e-8    1.58944e-7    2.57213e-7   4.05255e-7    6.21655e-7    9.28446e-7    1.35005e-6    1.9113e-6     2.63448e-6    3.53547e-6    4.6194e-6     5.87637e-6    7.27813e-6    8.77641e-6    1.03039e-5    1.17779e-5    1.31076e-5    1.42025e-5    1.49828e-5    1.53889e-5   1.53889e-5   1.49828e-5   1.42025e-5   1.31076e-5   1.17779e-5   1.03039e-5   8.77641e-6   7.27813e-6   5.87637e-6   4.6194e-6    3.53547e-6   2.63448e-6   1.9113e-6    1.35005e-6   9.28446e-7   6.21655e-7   4.05255e-7   2.57213e-7  1.58944e-7   9.56274e-8   5.60153e-8   3.1946e-8    1.77383e-8   9.58947e-9
  4.34312e-10   5.46896e-9    1.01163e-8    1.82191e-8    3.1946e-8     5.45372e-8    9.06473e-8    1.46691e-7   2.3112e-7     3.54535e-7    5.29501e-7    7.69945e-7    1.09003e-6    1.50247e-6    2.01631e-6    2.63448e-6    3.35135e-6    4.15078e-6    5.00526e-6    5.87637e-6    6.71705e-6    7.47539e-6    8.09982e-6    8.54482e-6    8.77641e-6   8.77641e-6   8.54482e-6   8.09982e-6   7.47539e-6   6.71705e-6   5.87637e-6   5.00526e-6   4.15078e-6   3.35135e-6   2.63448e-6   2.01631e-6   1.50247e-6   1.09003e-6   7.69945e-7   5.29501e-7   3.54535e-7   2.3112e-7    1.46691e-7  9.06473e-8   5.45372e-8   3.1946e-8    1.82191e-8   1.01163e-8   5.46896e-9
  4.54513e-10   3.03669e-9    5.61718e-9    1.01163e-8    1.77383e-8    3.02823e-8    5.03328e-8    8.14516e-8   1.28332e-7    1.96859e-7    2.9401e-7     4.2752e-7     6.05251e-7    8.3426e-7     1.11957e-6    1.46282e-6    1.86087e-6    2.30476e-6    2.77922e-6    3.26291e-6    3.72971e-6    4.15078e-6    4.4975e-6     4.74459e-6    4.87318e-6   4.87318e-6   4.74459e-6   4.4975e-6    4.15078e-6   3.72971e-6   3.26291e-6   2.77922e-6   2.30476e-6   1.86087e-6   1.46282e-6   1.11957e-6   8.3426e-7    6.05251e-7   4.2752e-7    2.9401e-7    1.96859e-7   1.28332e-7   8.14516e-8  5.03328e-8   3.02823e-8   1.77383e-8   1.01163e-8   5.61718e-9   3.03669e-9
  4.74713e-10   1.64166e-9    3.03669e-9    5.46896e-9    9.58947e-9    1.63708e-8    2.72103e-8    4.40334e-8   6.93772e-8    1.06424e-7    1.58944e-7    2.3112e-7     3.27203e-7    4.51007e-7    6.05251e-7    7.90813e-7    1.006e-6      1.24597e-6    1.50247e-6    1.76396e-6    2.01631e-6    2.24395e-6    2.43138e-6    2.56496e-6    2.63448e-6   2.63448e-6   2.56496e-6   2.43138e-6   2.24395e-6   2.01631e-6   1.76396e-6   1.50247e-6   1.24597e-6   1.006e-6     7.90813e-7   6.05251e-7   4.51007e-7   3.27203e-7   2.3112e-7    1.58944e-7   1.06424e-7   6.93772e-8   4.40334e-8  2.72103e-8   1.63708e-8   9.58947e-9   5.46896e-9   3.03669e-9   1.64166e-9

Now we can form our metadata we need to fully define our model. We will also fix the total flux to be the observed value 1.1. This is because total flux is degenerate with a global shift in the gain amplitudes making the problem degenerate. To fix this we use the observed total flux as our value.

julia
skymeta = (; ftot = 1.1, mimg = mimg ./ flux(mimg))
(ftot = 1.1, mimg = [1.6419398126916194e-9 3.0372105473444726e-9 5.4698920279304e-9 9.591108334634163e-9 1.6373632090899278e-8 2.7214944149646553e-8 4.404088682310447e-8 6.938903568369913e-8 1.0644170854928012e-7 1.5897140421510779e-7 2.3115983313045458e-7 3.272592438419348e-7 4.510842476340416e-7 6.053542398020764e-7 7.909475845251705e-7 1.0061714917927168e-6 1.2461849459797968e-6 1.502723858077153e-6 1.7642578193002676e-6 2.0166526840829715e-6 2.2443281409323587e-6 2.431799639052888e-6 2.5654019043369063e-6 2.6349308626683512e-6 2.6349308626683512e-6 2.5654019043369063e-6 2.431799639052888e-6 2.244328140932357e-6 2.0166526840829715e-6 1.7642578193002676e-6 1.502723858077153e-6 1.2461849459797968e-6 1.0061714917927168e-6 7.909475845251705e-7 6.053542398020754e-7 4.5108424763404237e-7 3.272592438419348e-7 2.3115983313045418e-7 1.5897140421510802e-7 1.0644170854928012e-7 6.938903568369913e-8 4.404088682310447e-8 2.72149441496466e-8 1.6373632090899278e-8 9.591108334634163e-9 5.4698920279304e-9 3.0372105473444726e-9 1.6419398126916194e-9; 3.0372105473444726e-9 5.618140103307805e-9 1.0118040644152322e-8 1.7741341777273437e-8 3.0287448845822705e-8 5.034137961561296e-8 8.146549888102353e-8 1.2835373709777678e-7 1.9689264940429196e-7 2.940604898402379e-7 4.275924597718894e-7 6.053542398020764e-7 8.344019823778814e-7 1.1197659425728733e-6 1.4630708918486186e-6 1.8611855585014145e-6 2.3051551784149337e-6 2.7796930899776574e-6 3.263470692162403e-6 3.7303430704848267e-6 4.151490236519409e-6 4.498269337078975e-6 4.745402761905639e-6 4.8740153236803395e-6 4.8740153236803395e-6 4.745402761905639e-6 4.498269337078975e-6 4.151490236519405e-6 3.7303430704848267e-6 3.263470692162403e-6 2.7796930899776574e-6 2.3051551784149337e-6 1.8611855585014145e-6 1.4630708918486186e-6 1.1197659425728733e-6 8.344019823778828e-7 6.053542398020764e-7 4.275924597718886e-7 2.940604898402384e-7 1.9689264940429196e-7 1.2835373709777678e-7 8.146549888102353e-8 5.034137961561296e-8 3.0287448845822705e-8 1.7741341777273437e-8 1.0118040644152322e-8 5.618140103307805e-9 3.0372105473444726e-9; 5.4698920279304e-9 1.0118040644152322e-8 1.822217755239726e-8 3.195143123585701e-8 5.4546457153905775e-8 9.066276662156884e-8 1.467160329962283e-7 2.3115983313045458e-7 3.5459561217323976e-7 5.295909203636775e-7 7.700765391304108e-7 1.0902182376728397e-6 1.502723858077153e-6 2.0166526840829736e-6 2.6349308626683538e-6 3.3519191015081937e-6 4.151490236519409e-6 5.00611361509195e-6 5.8773772987359084e-6 6.718195365326499e-6 7.476664193933783e-6 8.101199176950392e-6 8.546276371706097e-6 8.777902337498369e-6 8.777902337498369e-6 8.546276371706097e-6 8.101199176950392e-6 7.476664193933778e-6 6.718195365326499e-6 5.8773772987359084e-6 5.00611361509195e-6 4.151490236519409e-6 3.3519191015081937e-6 2.6349308626683538e-6 2.0166526840829715e-6 1.5027238580771545e-6 1.0902182376728397e-6 7.700765391304094e-7 5.295909203636784e-7 3.5459561217323976e-7 2.3115983313045458e-7 1.467160329962283e-7 9.066276662156917e-8 5.4546457153905775e-8 3.195143123585701e-8 1.822217755239726e-8 1.0118040644152322e-8 5.4698920279304e-9; 9.591108334634163e-9 1.7741341777273437e-8 3.195143123585701e-8 5.6024805766717734e-8 9.564374893731387e-8 1.5897140421510802e-7 2.572572474391242e-7 4.053240888209999e-7 6.217608892412182e-7 9.286040500087317e-7 1.350280311758629e-6 1.911628451260389e-6 2.6349308626683487e-6 3.53607242475841e-6 4.620183950447215e-6 5.8773772987359084e-6 7.279374511474409e-6 8.77790233749836e-6 1.0305607881811155e-5 1.1779928969907525e-5 1.3109855898349644e-5 1.420493833330359e-5 1.4985352932052472e-5 1.5391494354896698e-5 1.5391494354896698e-5 1.4985352932052472e-5 1.420493833330359e-5 1.3109855898349632e-5 1.1779928969907525e-5 1.0305607881811155e-5 8.77790233749836e-6 7.279374511474415e-6 5.8773772987359084e-6 4.620183950447215e-6 3.5360724247584035e-6 2.6349308626683512e-6 1.911628451260389e-6 1.3502803117586267e-6 9.286040500087334e-7 6.217608892412182e-7 4.053240888209999e-7 2.572572474391242e-7 1.5897140421510834e-7 9.564374893731387e-8 5.6024805766717734e-8 3.195143123585701e-8 1.7741341777273437e-8 9.591108334634163e-9; 1.6373632090899278e-8 3.0287448845822705e-8 5.4546457153905775e-8 9.564374893731387e-8 1.6327993619244742e-7 2.713908752539494e-7 4.391813099509043e-7 6.919562657809534e-7 1.0614502196981392e-6 1.5852830082272205e-6 2.3051551784149337e-6 3.263470692162409e-6 4.498269337078971e-6 6.0366692680024615e-6 7.887429643946967e-6 1.0033669748197654e-5 1.2427114358863923e-5 1.4985352932052472e-5 1.759340276874101e-5 2.011031637643879e-5 2.2380724912588598e-5 2.425021447245457e-5 2.5582513209206184e-5 2.627586480137954e-5 2.627586480137954e-5 2.5582513209206184e-5 2.425021447245457e-5 2.2380724912588574e-5 2.011031637643879e-5 1.759340276874101e-5 1.4985352932052472e-5 1.2427114358863933e-5 1.0033669748197654e-5 7.887429643946967e-6 6.0366692680024505e-6 4.498269337078975e-6 3.263470692162409e-6 2.305155178414932e-6 1.5852830082272232e-6 1.0614502196981392e-6 6.919562657809534e-7 4.391813099509043e-7 2.713908752539499e-7 1.6327993619244742e-7 9.564374893731387e-8 5.4546457153905775e-8 3.0287448845822705e-8 1.6373632090899278e-8; 2.7214944149646553e-8 5.034137961561296e-8 9.066276662156884e-8 1.5897140421510802e-7 2.713908752539494e-7 4.510842476340432e-7 7.299721134277668e-7 1.1501144659097052e-6 1.7642578193002719e-6 2.6349308626683538e-6 3.8314449163481485e-6 5.424280460691016e-6 7.476664193933783e-6 1.0033669748197674e-5 1.3109855898349667e-5 1.6677164864658932e-5 2.0655357427124005e-5 2.4907457358543318e-5 2.9242349595704554e-5 3.3425762468469194e-5 3.7199454289878896e-5 4.0306770594388894e-5 4.2521211196819014e-5 4.367364486286884e-5 4.367364486286884e-5 4.2521211196819014e-5 4.0306770594388894e-5 3.719945428987886e-5 3.3425762468469194e-5 2.9242349595704554e-5 2.4907457358543318e-5 2.065535742712402e-5 1.6677164864658932e-5 1.3109855898349667e-5 1.0033669748197654e-5 7.476664193933791e-6 5.424280460691016e-6 3.831444916348146e-6 2.634930862668358e-6 1.7642578193002719e-6 1.1501144659097052e-6 7.299721134277668e-7 4.510842476340439e-7 2.713908752539494e-7 1.5897140421510802e-7 9.066276662156884e-8 5.034137961561296e-8 2.7214944149646553e-8; 4.404088682310447e-8 8.146549888102353e-8 1.467160329962283e-7 2.572572474391242e-7 4.391813099509043e-7 7.299721134277668e-7 1.1812855119128452e-6 1.8611855585014124e-6 2.8550298879665212e-6 4.264006248603274e-6 6.200278457379807e-6 8.777902337498369e-6 1.2099195198373104e-5 1.6237097947765598e-5 2.1215170485438063e-5 2.6988007996490334e-5 3.342576246846914e-5 4.0306770594388826e-5 4.73217582922288e-5 5.409161281960394e-5 6.019843168713919e-5 6.522688094420994e-5 6.881042414061877e-5 7.067536277060657e-5 7.067536277060657e-5 6.881042414061877e-5 6.522688094420994e-5 6.0198431687139135e-5 5.409161281960394e-5 4.73217582922288e-5 4.0306770594388826e-5 3.342576246846914e-5 2.6988007996490334e-5 2.1215170485438063e-5 1.623709794776558e-5 1.2099195198373111e-5 8.777902337498369e-6 6.200278457379796e-6 4.264006248603281e-6 2.8550298879665212e-6 1.8611855585014124e-6 1.1812855119128452e-6 7.299721134277695e-7 4.391813099509043e-7 2.572572474391242e-7 1.467160329962283e-7 8.146549888102353e-8 4.404088682310447e-8; 6.938903568369913e-8 1.2835373709777678e-7 2.3115983313045458e-7 4.053240888209999e-7 6.919562657809534e-7 1.1501144659097052e-6 1.8611855585014124e-6 2.932408506022366e-6 4.498269337078975e-6 6.718195365326483e-6 9.76890735320735e-6 1.3830107031477007e-5 1.9063001404492862e-5 2.5582513209206164e-5 3.342576246846911e-5 4.252121119681889e-5 5.2664276130395184e-5 6.350571354983377e-5 7.455824375071895e-5 8.52245520669049e-5 9.484620827958101e-5 0.00010276882905545298 0.00010841491442445243 0.00011135323611774409 0.00011135323611774409 0.00010841491442445243 0.00010276882905545298 9.484620827958101e-5 8.52245520669049e-5 7.455824375071895e-5 6.350571354983377e-5 5.2664276130395184e-5 4.252121119681889e-5 3.342576246846911e-5 2.558251320920614e-5 1.906300140449289e-5 1.3830107031477007e-5 9.768907353207345e-6 6.7181953653264934e-6 4.498269337078975e-6 2.932408506022366e-6 1.8611855585014124e-6 1.150114465909707e-6 6.919562657809534e-7 4.053240888209999e-7 2.3115983313045458e-7 1.2835373709777678e-7 6.938903568369913e-8; 1.0644170854928012e-7 1.9689264940429196e-7 3.5459561217323976e-7 6.217608892412182e-7 1.0614502196981392e-6 1.7642578193002719e-6 2.8550298879665212e-6 4.498269337078975e-6 6.90027565646087e-6 1.0305607881811166e-5 1.49853529320525e-5 2.12151704854381e-5 2.924234959570453e-5 3.9243179965565206e-5 5.12746031365021e-5 6.522688094421001e-5 8.07861858231187e-5 9.741678330418785e-5 0.0001143711938502203 0.0001307331459929699 0.0001454926179507512 0.00015764579580172894 0.00016630680351526565 0.000170814143590246 0.000170814143590246 0.00016630680351526565 0.00015764579580172894 0.0001454926179507512 0.00013073314599296995 0.0001143711938502203 9.741678330418779e-5 8.07861858231187e-5 6.522688094421001e-5 5.12746031365021e-5 3.924317996556517e-5 2.924234959570458e-5 2.12151704854381e-5 1.4985352932052484e-5 1.0305607881811182e-5 6.90027565646087e-6 4.498269337078975e-6 2.8550298879665212e-6 1.7642578193002753e-6 1.0614502196981392e-6 6.217608892412182e-7 3.5459561217323976e-7 1.9689264940429196e-7 1.0644170854928012e-7; 1.5897140421510779e-7 2.940604898402379e-7 5.295909203636775e-7 9.286040500087317e-7 1.5852830082272205e-6 2.6349308626683538e-6 4.264006248603274e-6 6.718195365326483e-6 1.0305607881811166e-5 1.5391494354896675e-5 2.2380724912588574e-5 3.168499912956171e-5 4.3673644862868684e-5 5.860995196355548e-5 7.65789630049792e-5 9.741678330418758e-5 0.00012065470929131058 0.00014549261795075095 0.0001708141435902456 0.00019525082863865447 0.00021729419880417281 0.00023544505127523585 0.0002483803242702325 0.00025511206684333235 0.00025511206684333235 0.0002483803242702325 0.00023544505127523585 0.00021729419880417281 0.00019525082863865455 0.0001708141435902456 0.0001454926179507509 0.00012065470929131058 9.741678330418758e-5 7.65789630049792e-5 5.860995196355542e-5 4.367364486286876e-5 3.168499912956171e-5 2.2380724912588537e-5 1.5391494354896698e-5 1.0305607881811166e-5 6.718195365326483e-6 4.264006248603274e-6 2.634930862668358e-6 1.5852830082272205e-6 9.286040500087317e-7 5.295909203636775e-7 2.940604898402379e-7 1.5897140421510779e-7; 2.3115983313045458e-7 4.275924597718894e-7 7.700765391304108e-7 1.350280311758629e-6 2.3051551784149337e-6 3.8314449163481485e-6 6.200278457379807e-6 9.76890735320735e-6 1.49853529320525e-5 2.2380724912588574e-5 3.254374371086372e-5 4.607306042046038e-5 6.350571354983377e-5 8.5224552066905e-5 0.00011135323611774424 0.0001416534469446525 0.0001754436441188102 0.00021156037120802129 0.00024838032427023277 0.0002839136333325858 0.00031596683053652545 0.00034235992965458496 0.00036116907059275224 0.00037095767689938896 0.00037095767689938896 0.00036116907059275224 0.00034235992965458496 0.00031596683053652545 0.0002839136333325859 0.00024838032427023277 0.00021156037120802118 0.0001754436441188102 0.0001416534469446525 0.00011135323611774424 8.522455206690495e-5 6.350571354983388e-5 4.607306042046038e-5 3.2543743710863696e-5 2.238072491258862e-5 1.49853529320525e-5 9.76890735320735e-6 6.200278457379807e-6 3.831444916348156e-6 2.3051551784149337e-6 1.350280311758629e-6 7.700765391304108e-7 4.275924597718894e-7 2.3115983313045458e-7; 3.272592438419348e-7 6.053542398020764e-7 1.0902182376728397e-6 1.911628451260389e-6 3.263470692162409e-6 5.424280460691016e-6 8.777902337498369e-6 1.3830107031477007e-5 2.12151704854381e-5 3.168499912956171e-5 4.607306042046038e-5 6.522688094421008e-5 8.990676067944877e-5 0.0001206547092913108 0.00015764579580172894 0.0002005426258832297 0.0002483803242702331 0.0002995117541436506 0.00035163876009556546 0.0004019442292485252 0.00044732280967758225 0.00048468823576848475 0.0005113168466191734 0.0005251748419931728 0.0005251748419931728 0.0005113168466191734 0.00048468823576848475 0.00044732280967758225 0.0004019442292485254 0.00035163876009556546 0.0002995117541436505 0.0002483803242702331 0.0002005426258832297 0.00015764579580172894 0.00012065470929131066 8.990676067944894e-5 6.522688094421008e-5 4.60730604204603e-5 3.168499912956177e-5 2.12151704854381e-5 1.3830107031477007e-5 8.777902337498369e-6 5.424280460691025e-6 3.263470692162409e-6 1.911628451260389e-6 1.0902182376728397e-6 6.053542398020764e-7 3.272592438419348e-7; 4.510842476340416e-7 8.344019823778814e-7 1.502723858077153e-6 2.6349308626683487e-6 4.498269337078971e-6 7.476664193933783e-6 1.2099195198373104e-5 1.9063001404492862e-5 2.924234959570453e-5 4.3673644862868684e-5 6.350571354983377e-5 8.990676067944877e-5 0.0001239247607560054 0.00016630680351526535 0.0002172941988041729 0.00027642189248223184 0.0003423599296545847 0.0004128379467279341 0.00048468823576848383 0.0005540277735561485 0.0006165762368821491 0.0006680796105313319 0.0007047836826611451 0.0007238851245138184 0.0007238851245138184 0.0007047836826611451 0.0006680796105313319 0.0006165762368821491 0.0005540277735561485 0.00048468823576848383 0.0004128379467279341 0.0003423599296545847 0.00027642189248223184 0.0002172941988041729 0.00016630680351526516 0.00012392476075600558 8.990676067944877e-5 6.35057135498337e-5 4.367364486286876e-5 2.924234959570453e-5 1.9063001404492862e-5 1.2099195198373104e-5 7.476664193933798e-6 4.498269337078971e-6 2.6349308626683487e-6 1.502723858077153e-6 8.344019823778814e-7 4.510842476340416e-7; 6.053542398020764e-7 1.1197659425728733e-6 2.0166526840829736e-6 3.53607242475841e-6 6.0366692680024615e-6 1.0033669748197674e-5 1.6237097947765598e-5 2.5582513209206164e-5 3.9243179965565206e-5 5.860995196355548e-5 8.5224552066905e-5 0.0001206547092913108 0.00016630680351526535 0.00022318342780520376 0.00029160841953678227 0.00037095767689938896 0.00045944640284331504 0.0005540277735561492 0.0006504507307528019 0.0007435042643351424 0.0008274441883429026 0.0008965616220954052 0.0009458184201288619 0.0009714525203496029 0.0009714525203496029 0.0009458184201288619 0.0008965616220954052 0.0008274441883429026 0.0007435042643351424 0.0006504507307528019 0.0005540277735561492 0.00045944640284331504 0.00037095767689938896 0.00029160841953678227 0.0002231834278052034 0.0001663068035152656 0.0001206547092913108 8.522455206690495e-5 5.860995196355557e-5 3.9243179965565206e-5 2.5582513209206164e-5 1.6237097947765598e-5 1.003366974819769e-5 6.0366692680024615e-6 3.53607242475841e-6 2.0166526840829736e-6 1.1197659425728733e-6 6.053542398020764e-7; 7.909475845251705e-7 1.4630708918486186e-6 2.6349308626683538e-6 4.620183950447215e-6 7.887429643946967e-6 1.3109855898349667e-5 2.1215170485438063e-5 3.342576246846911e-5 5.12746031365021e-5 7.65789630049792e-5 0.00011135323611774424 0.00015764579580172894 0.0002172941988041729 0.00029160841953678227 0.0003810115794930779 0.0004846882357684846 0.0006003063969065673 0.0007238851245138195 0.0008498700438767397 0.0009714525203496029 0.0010811272789849299 0.001171435174231512 0.0012357934340146508 0.0012692865993680492 0.0012692865993680492 0.0012357934340146508 0.001171435174231512 0.0010811272789849299 0.0009714525203496029 0.0008498700438767397 0.0007238851245138195 0.0006003063969065675 0.0004846882357684846 0.0003810115794930779 0.00029160841953678184 0.0002172941988041732 0.00015764579580172894 0.00011135323611774417 7.657896300497933e-5 5.12746031365021e-5 3.342576246846911e-5 2.1215170485438063e-5 1.3109855898349693e-5 7.887429643946967e-6 4.620183950447215e-6 2.6349308626683538e-6 1.4630708918486186e-6 7.909475845251705e-7; 1.0061714917927168e-6 1.8611855585014145e-6 3.3519191015081937e-6 5.8773772987359084e-6 1.0033669748197654e-5 1.6677164864658932e-5 2.6988007996490334e-5 4.252121119681889e-5 6.522688094421001e-5 9.741678330418758e-5 0.0001416534469446525 0.0002005426258832297 0.00027642189248223184 0.00037095767689938896 0.0004846882357684846 0.0006165762368821498 0.0007636551330652295 0.0009208607370055721 0.0010811272789849292 0.0012357934340146508 0.0013753116747516273 0.0014901931554700534 0.001572063872976649 0.0016146708279858348 0.001614670827985835 0.001572063872976649 0.0014901931554700534 0.0013753116747516271 0.0012357934340146508 0.0010811272789849292 0.0009208607370055721 0.0007636551330652295 0.0006165762368821498 0.0004846882357684846 0.00037095767689938853 0.00027642189248223216 0.0002005426258832297 0.00014165344694465232 9.741678330418775e-5 6.522688094421001e-5 4.252121119681889e-5 2.6988007996490334e-5 1.667716486465896e-5 1.0033669748197654e-5 5.8773772987359084e-6 3.3519191015081937e-6 1.8611855585014145e-6 1.0061714917927168e-6; 1.2461849459797968e-6 2.3051551784149337e-6 4.151490236519409e-6 7.279374511474409e-6 1.2427114358863923e-5 2.0655357427124005e-5 3.342576246846914e-5 5.2664276130395184e-5 8.07861858231187e-5 0.00012065470929131058 0.0001754436441188102 0.0002483803242702331 0.0003423599296545847 0.00045944640284331504 0.0006003063969065673 0.0007636551330652295 0.0009458184201288615 0.001140524052967917 0.0013390207839805053 0.0015305812044682718 0.0017033803075180167 0.001845665765823035 0.001947066030594502 0.0019998365039776067 0.0019998365039776075 0.001947066030594502 0.001845665765823035 0.0017033803075180167 0.0015305812044682718 0.0013390207839805053 0.001140524052967917 0.0009458184201288619 0.0007636551330652295 0.0006003063969065673 0.00045944640284331455 0.0003423599296545851 0.0002483803242702331 0.00017544364411880986 0.0001206547092913108 8.07861858231187e-5 5.2664276130395184e-5 3.342576246846914e-5 2.065535742712404e-5 1.2427114358863923e-5 7.279374511474409e-6 4.151490236519409e-6 2.3051551784149337e-6 1.2461849459797968e-6; 1.502723858077153e-6 2.7796930899776574e-6 5.00611361509195e-6 8.77790233749836e-6 1.4985352932052472e-5 2.4907457358543318e-5 4.0306770594388826e-5 6.350571354983377e-5 9.741678330418785e-5 0.00014549261795075095 0.00021156037120802129 0.0002995117541436506 0.0004128379467279341 0.0005540277735561492 0.0007238851245138195 0.0009208607370055721 0.001140524052967917 0.001375311674751628 0.001614670827985835 0.0018456657658230359 0.002054037192267305 0.0022256134527109616 0.0023478879173311545 0.002411521689838708 0.002411521689838709 0.0023478879173311545 0.0022256134527109616 0.002054037192267305 0.0018456657658230363 0.001614670827985835 0.0013753116747516278 0.0011405240529679171 0.0009208607370055721 0.0007238851245138195 0.0005540277735561487 0.0004128379467279347 0.0002995117541436506 0.000211560371208021 0.0001454926179507512 9.741678330418785e-5 6.350571354983377e-5 4.0306770594388826e-5 2.4907457358543362e-5 1.4985352932052472e-5 8.77790233749836e-6 5.00611361509195e-6 2.7796930899776574e-6 1.502723858077153e-6; 1.7642578193002676e-6 3.263470692162403e-6 5.8773772987359084e-6 1.0305607881811155e-5 1.759340276874101e-5 2.9242349595704554e-5 4.73217582922288e-5 7.455824375071895e-5 0.0001143711938502203 0.0001708141435902456 0.00024838032427023277 0.00035163876009556546 0.00048468823576848383 0.0006504507307528019 0.0008498700438767397 0.0010811272789849292 0.0013390207839805053 0.001614670827985835 0.0018956880324740198 0.002166885314068744 0.0024115216898387077 0.0026129590713422886 0.0027565141757265126 0.0028312227658075133 0.0028312227658075133 0.0027565141757265126 0.0026129590713422886 0.0024115216898387077 0.002166885314068744 0.0018956880324740198 0.0016146708279858348 0.001339020783980506 0.0010811272789849292 0.0008498700438767397 0.0006504507307528012 0.0004846882357684846 0.00035163876009556546 0.00024838032427023244 0.0001708141435902459 0.0001143711938502203 7.455824375071895e-5 4.73217582922288e-5 2.9242349595704604e-5 1.759340276874101e-5 1.0305607881811155e-5 5.8773772987359084e-6 3.263470692162403e-6 1.7642578193002676e-6; 2.0166526840829715e-6 3.7303430704848267e-6 6.718195365326499e-6 1.1779928969907525e-5 2.011031637643879e-5 3.3425762468469194e-5 5.409161281960394e-5 8.52245520669049e-5 0.0001307331459929699 0.00019525082863865447 0.0002839136333325858 0.0004019442292485252 0.0005540277735561485 0.0007435042643351424 0.0009714525203496029 0.0012357934340146508 0.0015305812044682718 0.0018456657658230359 0.002166885314068744 0.002476880100466188 0.0027565141757265126 0.0029867692051444685 0.003150861314247364 0.0032362577212025905 0.0032362577212025905 0.003150861314247364 0.0029867692051444685 0.0027565141757265126 0.002476880100466188 0.002166885314068744 0.0018456657658230357 0.001530581204468272 0.0012357934340146508 0.0009714525203496029 0.0007435042643351416 0.0005540277735561492 0.0004019442292485252 0.0002839136333325856 0.00019525082863865482 0.0001307331459929699 8.52245520669049e-5 5.409161281960394e-5 3.3425762468469255e-5 2.011031637643879e-5 1.1779928969907525e-5 6.718195365326499e-6 3.7303430704848267e-6 2.0166526840829715e-6; 2.2443281409323587e-6 4.151490236519409e-6 7.476664193933783e-6 1.3109855898349644e-5 2.2380724912588598e-5 3.7199454289878896e-5 6.019843168713919e-5 9.484620827958101e-5 0.0001454926179507512 0.00021729419880417281 0.00031596683053652545 0.00044732280967758225 0.0006165762368821491 0.0008274441883429026 0.0010811272789849299 0.0013753116747516273 0.0017033803075180167 0.002054037192267305 0.0024115216898387077 0.0027565141757265126 0.0030677182959122987 0.0033239685893776383 0.0035065863207655477 0.0036016237859556825 0.0036016237859556825 0.0035065863207655477 0.0033239685893776383 0.0030677182959122987 0.0027565141757265135 0.0024115216898387077 0.002054037192267305 0.001703380307518017 0.0013753116747516273 0.0010811272789849299 0.0008274441883429017 0.0006165762368821499 0.00044732280967758225 0.00031596683053652507 0.0002172941988041732 0.0001454926179507512 9.484620827958101e-5 6.019843168713919e-5 3.719945428987896e-5 2.2380724912588598e-5 1.3109855898349644e-5 7.476664193933783e-6 4.151490236519409e-6 2.2443281409323587e-6; 2.431799639052888e-6 4.498269337078975e-6 8.101199176950392e-6 1.420493833330359e-5 2.425021447245457e-5 4.0306770594388894e-5 6.522688094420994e-5 0.00010276882905545298 0.00015764579580172894 0.00023544505127523585 0.00034235992965458496 0.00048468823576848475 0.0006680796105313319 0.0008965616220954052 0.001171435174231512 0.0014901931554700534 0.001845665765823035 0.0022256134527109616 0.0026129590713422886 0.0029867692051444685 0.0033239685893776383 0.003601623785955681 0.003799495801715948 0.0039024718636076175 0.003902471863607618 0.003799495801715948 0.003601623785955681 0.0033239685893776375 0.0029867692051444685 0.0026129590713422886 0.002225613452710961 0.0018456657658230359 0.0014901931554700534 0.001171435174231512 0.0008965616220954041 0.0006680796105313329 0.00048468823576848475 0.0003423599296545846 0.0002354450512752363 0.00015764579580172894 0.00010276882905545298 6.522688094420994e-5 4.0306770594388975e-5 2.425021447245457e-5 1.420493833330359e-5 8.101199176950392e-6 4.498269337078975e-6 2.431799639052888e-6; 2.5654019043369063e-6 4.745402761905639e-6 8.546276371706097e-6 1.4985352932052472e-5 2.5582513209206184e-5 4.2521211196819014e-5 6.881042414061877e-5 0.00010841491442445243 0.00016630680351526565 0.0002483803242702325 0.00036116907059275224 0.0005113168466191734 0.0007047836826611451 0.0009458184201288619 0.0012357934340146508 0.001572063872976649 0.001947066030594502 0.0023478879173311545 0.0027565141757265126 0.003150861314247364 0.0035065863207655477 0.003799495801715948 0.0040082388403669756 0.0041168723729309055 0.0041168723729309055 0.0040082388403669756 0.003799495801715948 0.003506586320765547 0.003150861314247365 0.0027565141757265126 0.002347887917331154 0.001947066030594503 0.001572063872976649 0.0012357934340146508 0.0009458184201288609 0.0007047836826611461 0.0005113168466191734 0.0003611690705927517 0.00024838032427023293 0.00016630680351526565 0.00010841491442445243 6.881042414061877e-5 4.252121119681908e-5 2.5582513209206184e-5 1.4985352932052472e-5 8.546276371706097e-6 4.745402761905639e-6 2.5654019043369063e-6; 2.6349308626683512e-6 4.8740153236803395e-6 8.777902337498369e-6 1.5391494354896698e-5 2.627586480137954e-5 4.367364486286884e-5 7.067536277060657e-5 0.00011135323611774409 0.000170814143590246 0.00025511206684333235 0.00037095767689938896 0.0005251748419931728 0.0007238851245138184 0.0009714525203496029 0.0012692865993680492 0.0016146708279858348 0.0019998365039776067 0.002411521689838708 0.0028312227658075133 0.0032362577212025905 0.0036016237859556825 0.0039024718636076175 0.0041168723729309055 0.00422845015229931 0.00422845015229931 0.0041168723729309055 0.0039024718636076175 0.003601623785955681 0.0032362577212025905 0.0028312227658075133 0.0024115216898387077 0.0019998365039776075 0.0016146708279858348 0.0012692865993680492 0.0009714525203496018 0.0007238851245138195 0.0005251748419931728 0.00037095767689938853 0.00025511206684333273 0.000170814143590246 0.00011135323611774409 7.067536277060657e-5 4.367364486286892e-5 2.627586480137954e-5 1.5391494354896698e-5 8.777902337498369e-6 4.8740153236803395e-6 2.6349308626683512e-6; 2.6349308626683512e-6 4.8740153236803395e-6 8.777902337498369e-6 1.5391494354896698e-5 2.627586480137954e-5 4.367364486286884e-5 7.067536277060657e-5 0.00011135323611774409 0.000170814143590246 0.00025511206684333235 0.00037095767689938896 0.0005251748419931728 0.0007238851245138184 0.0009714525203496029 0.0012692865993680492 0.001614670827985835 0.0019998365039776075 0.002411521689838709 0.0028312227658075133 0.0032362577212025905 0.0036016237859556825 0.003902471863607618 0.0041168723729309055 0.00422845015229931 0.00422845015229931 0.0041168723729309055 0.003902471863607618 0.003601623785955681 0.0032362577212025905 0.0028312227658075133 0.002411521689838708 0.0019998365039776075 0.001614670827985835 0.0012692865993680492 0.0009714525203496018 0.0007238851245138195 0.0005251748419931728 0.00037095767689938853 0.00025511206684333273 0.000170814143590246 0.00011135323611774409 7.067536277060657e-5 4.367364486286892e-5 2.627586480137954e-5 1.5391494354896698e-5 8.777902337498369e-6 4.8740153236803395e-6 2.6349308626683512e-6; 2.5654019043369063e-6 4.745402761905639e-6 8.546276371706097e-6 1.4985352932052472e-5 2.5582513209206184e-5 4.2521211196819014e-5 6.881042414061877e-5 0.00010841491442445243 0.00016630680351526565 0.0002483803242702325 0.00036116907059275224 0.0005113168466191734 0.0007047836826611451 0.0009458184201288619 0.0012357934340146508 0.001572063872976649 0.001947066030594502 0.0023478879173311545 0.0027565141757265126 0.003150861314247364 0.0035065863207655477 0.003799495801715948 0.0040082388403669756 0.0041168723729309055 0.0041168723729309055 0.0040082388403669756 0.003799495801715948 0.003506586320765547 0.003150861314247365 0.0027565141757265126 0.002347887917331154 0.001947066030594503 0.001572063872976649 0.0012357934340146508 0.0009458184201288609 0.0007047836826611461 0.0005113168466191734 0.0003611690705927517 0.00024838032427023293 0.00016630680351526565 0.00010841491442445243 6.881042414061877e-5 4.252121119681908e-5 2.5582513209206184e-5 1.4985352932052472e-5 8.546276371706097e-6 4.745402761905639e-6 2.5654019043369063e-6; 2.431799639052888e-6 4.498269337078975e-6 8.101199176950392e-6 1.420493833330359e-5 2.425021447245457e-5 4.0306770594388894e-5 6.522688094420994e-5 0.00010276882905545298 0.00015764579580172894 0.00023544505127523585 0.00034235992965458496 0.00048468823576848475 0.0006680796105313319 0.0008965616220954052 0.001171435174231512 0.0014901931554700534 0.001845665765823035 0.0022256134527109616 0.0026129590713422886 0.0029867692051444685 0.0033239685893776383 0.003601623785955681 0.003799495801715948 0.0039024718636076175 0.003902471863607618 0.003799495801715948 0.003601623785955681 0.0033239685893776375 0.0029867692051444685 0.0026129590713422886 0.002225613452710961 0.0018456657658230359 0.0014901931554700534 0.001171435174231512 0.0008965616220954041 0.0006680796105313329 0.00048468823576848475 0.0003423599296545846 0.0002354450512752363 0.00015764579580172894 0.00010276882905545298 6.522688094420994e-5 4.0306770594388975e-5 2.425021447245457e-5 1.420493833330359e-5 8.101199176950392e-6 4.498269337078975e-6 2.431799639052888e-6; 2.244328140932357e-6 4.151490236519405e-6 7.476664193933778e-6 1.3109855898349632e-5 2.2380724912588574e-5 3.719945428987886e-5 6.0198431687139135e-5 9.484620827958101e-5 0.0001454926179507512 0.00021729419880417281 0.00031596683053652545 0.00044732280967758225 0.0006165762368821491 0.0008274441883429026 0.0010811272789849299 0.0013753116747516271 0.0017033803075180167 0.002054037192267305 0.0024115216898387077 0.0027565141757265126 0.0030677182959122987 0.0033239685893776375 0.003506586320765547 0.003601623785955681 0.003601623785955681 0.003506586320765547 0.0033239685893776375 0.0030677182959122974 0.0027565141757265126 0.0024115216898387077 0.0020540371922673044 0.001703380307518017 0.0013753116747516271 0.0010811272789849299 0.0008274441883429017 0.0006165762368821499 0.00044732280967758225 0.00031596683053652507 0.0002172941988041732 0.0001454926179507512 9.484620827958101e-5 6.0198431687139135e-5 3.7199454289878936e-5 2.2380724912588574e-5 1.3109855898349632e-5 7.476664193933778e-6 4.151490236519405e-6 2.244328140932357e-6; 2.0166526840829715e-6 3.7303430704848267e-6 6.718195365326499e-6 1.1779928969907525e-5 2.011031637643879e-5 3.3425762468469194e-5 5.409161281960394e-5 8.52245520669049e-5 0.00013073314599296995 0.00019525082863865455 0.0002839136333325859 0.0004019442292485254 0.0005540277735561485 0.0007435042643351424 0.0009714525203496029 0.0012357934340146508 0.0015305812044682718 0.0018456657658230363 0.002166885314068744 0.002476880100466188 0.0027565141757265135 0.0029867692051444685 0.003150861314247365 0.0032362577212025905 0.0032362577212025905 0.003150861314247365 0.0029867692051444685 0.0027565141757265126 0.002476880100466188 0.002166885314068744 0.0018456657658230359 0.0015305812044682727 0.0012357934340146508 0.0009714525203496029 0.0007435042643351416 0.0005540277735561492 0.0004019442292485254 0.0002839136333325856 0.00019525082863865493 0.00013073314599296995 8.52245520669049e-5 5.409161281960394e-5 3.3425762468469255e-5 2.011031637643879e-5 1.1779928969907525e-5 6.718195365326499e-6 3.7303430704848267e-6 2.0166526840829715e-6; 1.7642578193002676e-6 3.263470692162403e-6 5.8773772987359084e-6 1.0305607881811155e-5 1.759340276874101e-5 2.9242349595704554e-5 4.73217582922288e-5 7.455824375071895e-5 0.0001143711938502203 0.0001708141435902456 0.00024838032427023277 0.00035163876009556546 0.00048468823576848383 0.0006504507307528019 0.0008498700438767397 0.0010811272789849292 0.0013390207839805053 0.001614670827985835 0.0018956880324740198 0.002166885314068744 0.0024115216898387077 0.0026129590713422886 0.0027565141757265126 0.0028312227658075133 0.0028312227658075133 0.0027565141757265126 0.0026129590713422886 0.0024115216898387077 0.002166885314068744 0.0018956880324740198 0.0016146708279858348 0.001339020783980506 0.0010811272789849292 0.0008498700438767397 0.0006504507307528012 0.0004846882357684846 0.00035163876009556546 0.00024838032427023244 0.0001708141435902459 0.0001143711938502203 7.455824375071895e-5 4.73217582922288e-5 2.9242349595704604e-5 1.759340276874101e-5 1.0305607881811155e-5 5.8773772987359084e-6 3.263470692162403e-6 1.7642578193002676e-6; 1.502723858077153e-6 2.7796930899776574e-6 5.00611361509195e-6 8.77790233749836e-6 1.4985352932052472e-5 2.4907457358543318e-5 4.0306770594388826e-5 6.350571354983377e-5 9.741678330418779e-5 0.0001454926179507509 0.00021156037120802118 0.0002995117541436505 0.0004128379467279341 0.0005540277735561492 0.0007238851245138195 0.0009208607370055721 0.001140524052967917 0.0013753116747516278 0.0016146708279858348 0.0018456657658230357 0.002054037192267305 0.002225613452710961 0.002347887917331154 0.0024115216898387077 0.002411521689838708 0.002347887917331154 0.002225613452710961 0.0020540371922673044 0.0018456657658230359 0.0016146708279858348 0.0013753116747516273 0.001140524052967917 0.0009208607370055721 0.0007238851245138195 0.0005540277735561485 0.0004128379467279347 0.0002995117541436505 0.0002115603712080209 0.00014549261795075114 9.741678330418779e-5 6.350571354983377e-5 4.0306770594388826e-5 2.4907457358543362e-5 1.4985352932052472e-5 8.77790233749836e-6 5.00611361509195e-6 2.7796930899776574e-6 1.502723858077153e-6; 1.2461849459797968e-6 2.3051551784149337e-6 4.151490236519409e-6 7.279374511474415e-6 1.2427114358863933e-5 2.065535742712402e-5 3.342576246846914e-5 5.2664276130395184e-5 8.07861858231187e-5 0.00012065470929131058 0.0001754436441188102 0.0002483803242702331 0.0003423599296545847 0.00045944640284331504 0.0006003063969065675 0.0007636551330652295 0.0009458184201288619 0.0011405240529679171 0.001339020783980506 0.001530581204468272 0.001703380307518017 0.0018456657658230359 0.001947066030594503 0.0019998365039776075 0.0019998365039776075 0.001947066030594503 0.0018456657658230359 0.001703380307518017 0.0015305812044682727 0.001339020783980506 0.001140524052967917 0.0009458184201288619 0.0007636551330652295 0.0006003063969065675 0.00045944640284331455 0.0003423599296545853 0.0002483803242702331 0.00017544364411880995 0.0001206547092913108 8.07861858231187e-5 5.2664276130395184e-5 3.342576246846914e-5 2.065535742712406e-5 1.2427114358863933e-5 7.279374511474415e-6 4.151490236519409e-6 2.3051551784149337e-6 1.2461849459797968e-6; 1.0061714917927168e-6 1.8611855585014145e-6 3.3519191015081937e-6 5.8773772987359084e-6 1.0033669748197654e-5 1.6677164864658932e-5 2.6988007996490334e-5 4.252121119681889e-5 6.522688094421001e-5 9.741678330418758e-5 0.0001416534469446525 0.0002005426258832297 0.00027642189248223184 0.00037095767689938896 0.0004846882357684846 0.0006165762368821498 0.0007636551330652295 0.0009208607370055721 0.0010811272789849292 0.0012357934340146508 0.0013753116747516273 0.0014901931554700534 0.001572063872976649 0.0016146708279858348 0.001614670827985835 0.001572063872976649 0.0014901931554700534 0.0013753116747516271 0.0012357934340146508 0.0010811272789849292 0.0009208607370055721 0.0007636551330652295 0.0006165762368821498 0.0004846882357684846 0.00037095767689938853 0.00027642189248223216 0.0002005426258832297 0.00014165344694465232 9.741678330418775e-5 6.522688094421001e-5 4.252121119681889e-5 2.6988007996490334e-5 1.667716486465896e-5 1.0033669748197654e-5 5.8773772987359084e-6 3.3519191015081937e-6 1.8611855585014145e-6 1.0061714917927168e-6; 7.909475845251705e-7 1.4630708918486186e-6 2.6349308626683538e-6 4.620183950447215e-6 7.887429643946967e-6 1.3109855898349667e-5 2.1215170485438063e-5 3.342576246846911e-5 5.12746031365021e-5 7.65789630049792e-5 0.00011135323611774424 0.00015764579580172894 0.0002172941988041729 0.00029160841953678227 0.0003810115794930779 0.0004846882357684846 0.0006003063969065673 0.0007238851245138195 0.0008498700438767397 0.0009714525203496029 0.0010811272789849299 0.001171435174231512 0.0012357934340146508 0.0012692865993680492 0.0012692865993680492 0.0012357934340146508 0.001171435174231512 0.0010811272789849299 0.0009714525203496029 0.0008498700438767397 0.0007238851245138195 0.0006003063969065675 0.0004846882357684846 0.0003810115794930779 0.00029160841953678184 0.0002172941988041732 0.00015764579580172894 0.00011135323611774417 7.657896300497933e-5 5.12746031365021e-5 3.342576246846911e-5 2.1215170485438063e-5 1.3109855898349693e-5 7.887429643946967e-6 4.620183950447215e-6 2.6349308626683538e-6 1.4630708918486186e-6 7.909475845251705e-7; 6.053542398020754e-7 1.1197659425728733e-6 2.0166526840829715e-6 3.5360724247584035e-6 6.0366692680024505e-6 1.0033669748197654e-5 1.623709794776558e-5 2.558251320920614e-5 3.924317996556517e-5 5.860995196355542e-5 8.522455206690495e-5 0.00012065470929131066 0.00016630680351526516 0.0002231834278052034 0.00029160841953678184 0.00037095767689938853 0.00045944640284331455 0.0005540277735561487 0.0006504507307528012 0.0007435042643351416 0.0008274441883429017 0.0008965616220954041 0.0009458184201288609 0.0009714525203496018 0.0009714525203496018 0.0009458184201288609 0.0008965616220954041 0.0008274441883429017 0.0007435042643351416 0.0006504507307528012 0.0005540277735561485 0.00045944640284331455 0.00037095767689938853 0.00029160841953678184 0.00022318342780520322 0.00016630680351526544 0.00012065470929131066 8.522455206690484e-5 5.8609951963555525e-5 3.924317996556517e-5 2.558251320920614e-5 1.623709794776558e-5 1.0033669748197674e-5 6.0366692680024505e-6 3.5360724247584035e-6 2.0166526840829715e-6 1.1197659425728733e-6 6.053542398020754e-7; 4.5108424763404237e-7 8.344019823778828e-7 1.5027238580771545e-6 2.6349308626683512e-6 4.498269337078975e-6 7.476664193933791e-6 1.2099195198373111e-5 1.906300140449289e-5 2.924234959570458e-5 4.367364486286876e-5 6.350571354983388e-5 8.990676067944894e-5 0.00012392476075600558 0.0001663068035152656 0.0002172941988041732 0.00027642189248223216 0.0003423599296545851 0.0004128379467279347 0.0004846882357684846 0.0005540277735561492 0.0006165762368821499 0.0006680796105313329 0.0007047836826611461 0.0007238851245138195 0.0007238851245138195 0.0007047836826611461 0.0006680796105313329 0.0006165762368821499 0.0005540277735561492 0.0004846882357684846 0.0004128379467279347 0.0003423599296545853 0.00027642189248223216 0.0002172941988041732 0.00016630680351526544 0.00012392476075600575 8.990676067944894e-5 6.350571354983377e-5 4.367364486286884e-5 2.924234959570458e-5 1.906300140449289e-5 1.2099195198373111e-5 7.476664193933804e-6 4.498269337078975e-6 2.6349308626683512e-6 1.5027238580771545e-6 8.344019823778828e-7 4.5108424763404237e-7; 3.272592438419348e-7 6.053542398020764e-7 1.0902182376728397e-6 1.911628451260389e-6 3.263470692162409e-6 5.424280460691016e-6 8.777902337498369e-6 1.3830107031477007e-5 2.12151704854381e-5 3.168499912956171e-5 4.607306042046038e-5 6.522688094421008e-5 8.990676067944877e-5 0.0001206547092913108 0.00015764579580172894 0.0002005426258832297 0.0002483803242702331 0.0002995117541436506 0.00035163876009556546 0.0004019442292485252 0.00044732280967758225 0.00048468823576848475 0.0005113168466191734 0.0005251748419931728 0.0005251748419931728 0.0005113168466191734 0.00048468823576848475 0.00044732280967758225 0.0004019442292485254 0.00035163876009556546 0.0002995117541436505 0.0002483803242702331 0.0002005426258832297 0.00015764579580172894 0.00012065470929131066 8.990676067944894e-5 6.522688094421008e-5 4.60730604204603e-5 3.168499912956177e-5 2.12151704854381e-5 1.3830107031477007e-5 8.777902337498369e-6 5.424280460691025e-6 3.263470692162409e-6 1.911628451260389e-6 1.0902182376728397e-6 6.053542398020764e-7 3.272592438419348e-7; 2.3115983313045418e-7 4.275924597718886e-7 7.700765391304094e-7 1.3502803117586267e-6 2.305155178414932e-6 3.831444916348146e-6 6.200278457379796e-6 9.768907353207345e-6 1.4985352932052484e-5 2.2380724912588537e-5 3.2543743710863696e-5 4.60730604204603e-5 6.35057135498337e-5 8.522455206690495e-5 0.00011135323611774417 0.00014165344694465232 0.00017544364411880986 0.000211560371208021 0.00024838032427023244 0.0002839136333325856 0.00031596683053652507 0.0003423599296545846 0.0003611690705927517 0.00037095767689938853 0.00037095767689938853 0.0003611690705927517 0.0003423599296545846 0.00031596683053652507 0.0002839136333325856 0.00024838032427023244 0.0002115603712080209 0.00017544364411880995 0.00014165344694465232 0.00011135323611774417 8.522455206690484e-5 6.350571354983377e-5 4.60730604204603e-5 3.2543743710863635e-5 2.2380724912588574e-5 1.4985352932052484e-5 9.768907353207345e-6 6.200278457379796e-6 3.831444916348152e-6 2.305155178414932e-6 1.3502803117586267e-6 7.700765391304094e-7 4.275924597718886e-7 2.3115983313045418e-7; 1.5897140421510802e-7 2.940604898402384e-7 5.295909203636784e-7 9.286040500087334e-7 1.5852830082272232e-6 2.634930862668358e-6 4.264006248603281e-6 6.7181953653264934e-6 1.0305607881811182e-5 1.5391494354896698e-5 2.238072491258862e-5 3.168499912956177e-5 4.367364486286876e-5 5.860995196355557e-5 7.657896300497933e-5 9.741678330418775e-5 0.0001206547092913108 0.0001454926179507512 0.0001708141435902459 0.00019525082863865482 0.0002172941988041732 0.0002354450512752363 0.00024838032427023293 0.00025511206684333273 0.00025511206684333273 0.00024838032427023293 0.0002354450512752363 0.0002172941988041732 0.00019525082863865493 0.0001708141435902459 0.00014549261795075114 0.0001206547092913108 9.741678330418775e-5 7.657896300497933e-5 5.8609951963555525e-5 4.367364486286884e-5 3.168499912956177e-5 2.2380724912588574e-5 1.539149435489673e-5 1.0305607881811182e-5 6.7181953653264934e-6 4.264006248603281e-6 2.634930862668363e-6 1.5852830082272232e-6 9.286040500087334e-7 5.295909203636784e-7 2.940604898402384e-7 1.5897140421510802e-7; 1.0644170854928012e-7 1.9689264940429196e-7 3.5459561217323976e-7 6.217608892412182e-7 1.0614502196981392e-6 1.7642578193002719e-6 2.8550298879665212e-6 4.498269337078975e-6 6.90027565646087e-6 1.0305607881811166e-5 1.49853529320525e-5 2.12151704854381e-5 2.924234959570453e-5 3.9243179965565206e-5 5.12746031365021e-5 6.522688094421001e-5 8.07861858231187e-5 9.741678330418785e-5 0.0001143711938502203 0.0001307331459929699 0.0001454926179507512 0.00015764579580172894 0.00016630680351526565 0.000170814143590246 0.000170814143590246 0.00016630680351526565 0.00015764579580172894 0.0001454926179507512 0.00013073314599296995 0.0001143711938502203 9.741678330418779e-5 8.07861858231187e-5 6.522688094421001e-5 5.12746031365021e-5 3.924317996556517e-5 2.924234959570458e-5 2.12151704854381e-5 1.4985352932052484e-5 1.0305607881811182e-5 6.90027565646087e-6 4.498269337078975e-6 2.8550298879665212e-6 1.7642578193002753e-6 1.0614502196981392e-6 6.217608892412182e-7 3.5459561217323976e-7 1.9689264940429196e-7 1.0644170854928012e-7; 6.938903568369913e-8 1.2835373709777678e-7 2.3115983313045458e-7 4.053240888209999e-7 6.919562657809534e-7 1.1501144659097052e-6 1.8611855585014124e-6 2.932408506022366e-6 4.498269337078975e-6 6.718195365326483e-6 9.76890735320735e-6 1.3830107031477007e-5 1.9063001404492862e-5 2.5582513209206164e-5 3.342576246846911e-5 4.252121119681889e-5 5.2664276130395184e-5 6.350571354983377e-5 7.455824375071895e-5 8.52245520669049e-5 9.484620827958101e-5 0.00010276882905545298 0.00010841491442445243 0.00011135323611774409 0.00011135323611774409 0.00010841491442445243 0.00010276882905545298 9.484620827958101e-5 8.52245520669049e-5 7.455824375071895e-5 6.350571354983377e-5 5.2664276130395184e-5 4.252121119681889e-5 3.342576246846911e-5 2.558251320920614e-5 1.906300140449289e-5 1.3830107031477007e-5 9.768907353207345e-6 6.7181953653264934e-6 4.498269337078975e-6 2.932408506022366e-6 1.8611855585014124e-6 1.150114465909707e-6 6.919562657809534e-7 4.053240888209999e-7 2.3115983313045458e-7 1.2835373709777678e-7 6.938903568369913e-8; 4.404088682310447e-8 8.146549888102353e-8 1.467160329962283e-7 2.572572474391242e-7 4.391813099509043e-7 7.299721134277668e-7 1.1812855119128452e-6 1.8611855585014124e-6 2.8550298879665212e-6 4.264006248603274e-6 6.200278457379807e-6 8.777902337498369e-6 1.2099195198373104e-5 1.6237097947765598e-5 2.1215170485438063e-5 2.6988007996490334e-5 3.342576246846914e-5 4.0306770594388826e-5 4.73217582922288e-5 5.409161281960394e-5 6.019843168713919e-5 6.522688094420994e-5 6.881042414061877e-5 7.067536277060657e-5 7.067536277060657e-5 6.881042414061877e-5 6.522688094420994e-5 6.0198431687139135e-5 5.409161281960394e-5 4.73217582922288e-5 4.0306770594388826e-5 3.342576246846914e-5 2.6988007996490334e-5 2.1215170485438063e-5 1.623709794776558e-5 1.2099195198373111e-5 8.777902337498369e-6 6.200278457379796e-6 4.264006248603281e-6 2.8550298879665212e-6 1.8611855585014124e-6 1.1812855119128452e-6 7.299721134277695e-7 4.391813099509043e-7 2.572572474391242e-7 1.467160329962283e-7 8.146549888102353e-8 4.404088682310447e-8; 2.72149441496466e-8 5.034137961561296e-8 9.066276662156917e-8 1.5897140421510834e-7 2.713908752539499e-7 4.510842476340439e-7 7.299721134277695e-7 1.150114465909707e-6 1.7642578193002753e-6 2.634930862668358e-6 3.831444916348156e-6 5.424280460691025e-6 7.476664193933798e-6 1.003366974819769e-5 1.3109855898349693e-5 1.667716486465896e-5 2.065535742712404e-5 2.4907457358543362e-5 2.9242349595704604e-5 3.3425762468469255e-5 3.719945428987896e-5 4.0306770594388975e-5 4.252121119681908e-5 4.367364486286892e-5 4.367364486286892e-5 4.252121119681908e-5 4.0306770594388975e-5 3.7199454289878936e-5 3.3425762468469255e-5 2.9242349595704604e-5 2.4907457358543362e-5 2.065535742712406e-5 1.667716486465896e-5 1.3109855898349693e-5 1.0033669748197674e-5 7.476664193933804e-6 5.424280460691025e-6 3.831444916348152e-6 2.634930862668363e-6 1.7642578193002753e-6 1.150114465909707e-6 7.299721134277695e-7 4.5108424763404475e-7 2.713908752539499e-7 1.5897140421510834e-7 9.066276662156917e-8 5.034137961561296e-8 2.72149441496466e-8; 1.6373632090899278e-8 3.0287448845822705e-8 5.4546457153905775e-8 9.564374893731387e-8 1.6327993619244742e-7 2.713908752539494e-7 4.391813099509043e-7 6.919562657809534e-7 1.0614502196981392e-6 1.5852830082272205e-6 2.3051551784149337e-6 3.263470692162409e-6 4.498269337078971e-6 6.0366692680024615e-6 7.887429643946967e-6 1.0033669748197654e-5 1.2427114358863923e-5 1.4985352932052472e-5 1.759340276874101e-5 2.011031637643879e-5 2.2380724912588598e-5 2.425021447245457e-5 2.5582513209206184e-5 2.627586480137954e-5 2.627586480137954e-5 2.5582513209206184e-5 2.425021447245457e-5 2.2380724912588574e-5 2.011031637643879e-5 1.759340276874101e-5 1.4985352932052472e-5 1.2427114358863933e-5 1.0033669748197654e-5 7.887429643946967e-6 6.0366692680024505e-6 4.498269337078975e-6 3.263470692162409e-6 2.305155178414932e-6 1.5852830082272232e-6 1.0614502196981392e-6 6.919562657809534e-7 4.391813099509043e-7 2.713908752539499e-7 1.6327993619244742e-7 9.564374893731387e-8 5.4546457153905775e-8 3.0287448845822705e-8 1.6373632090899278e-8; 9.591108334634163e-9 1.7741341777273437e-8 3.195143123585701e-8 5.6024805766717734e-8 9.564374893731387e-8 1.5897140421510802e-7 2.572572474391242e-7 4.053240888209999e-7 6.217608892412182e-7 9.286040500087317e-7 1.350280311758629e-6 1.911628451260389e-6 2.6349308626683487e-6 3.53607242475841e-6 4.620183950447215e-6 5.8773772987359084e-6 7.279374511474409e-6 8.77790233749836e-6 1.0305607881811155e-5 1.1779928969907525e-5 1.3109855898349644e-5 1.420493833330359e-5 1.4985352932052472e-5 1.5391494354896698e-5 1.5391494354896698e-5 1.4985352932052472e-5 1.420493833330359e-5 1.3109855898349632e-5 1.1779928969907525e-5 1.0305607881811155e-5 8.77790233749836e-6 7.279374511474415e-6 5.8773772987359084e-6 4.620183950447215e-6 3.5360724247584035e-6 2.6349308626683512e-6 1.911628451260389e-6 1.3502803117586267e-6 9.286040500087334e-7 6.217608892412182e-7 4.053240888209999e-7 2.572572474391242e-7 1.5897140421510834e-7 9.564374893731387e-8 5.6024805766717734e-8 3.195143123585701e-8 1.7741341777273437e-8 9.591108334634163e-9; 5.4698920279304e-9 1.0118040644152322e-8 1.822217755239726e-8 3.195143123585701e-8 5.4546457153905775e-8 9.066276662156884e-8 1.467160329962283e-7 2.3115983313045458e-7 3.5459561217323976e-7 5.295909203636775e-7 7.700765391304108e-7 1.0902182376728397e-6 1.502723858077153e-6 2.0166526840829736e-6 2.6349308626683538e-6 3.3519191015081937e-6 4.151490236519409e-6 5.00611361509195e-6 5.8773772987359084e-6 6.718195365326499e-6 7.476664193933783e-6 8.101199176950392e-6 8.546276371706097e-6 8.777902337498369e-6 8.777902337498369e-6 8.546276371706097e-6 8.101199176950392e-6 7.476664193933778e-6 6.718195365326499e-6 5.8773772987359084e-6 5.00611361509195e-6 4.151490236519409e-6 3.3519191015081937e-6 2.6349308626683538e-6 2.0166526840829715e-6 1.5027238580771545e-6 1.0902182376728397e-6 7.700765391304094e-7 5.295909203636784e-7 3.5459561217323976e-7 2.3115983313045458e-7 1.467160329962283e-7 9.066276662156917e-8 5.4546457153905775e-8 3.195143123585701e-8 1.822217755239726e-8 1.0118040644152322e-8 5.4698920279304e-9; 3.0372105473444726e-9 5.618140103307805e-9 1.0118040644152322e-8 1.7741341777273437e-8 3.0287448845822705e-8 5.034137961561296e-8 8.146549888102353e-8 1.2835373709777678e-7 1.9689264940429196e-7 2.940604898402379e-7 4.275924597718894e-7 6.053542398020764e-7 8.344019823778814e-7 1.1197659425728733e-6 1.4630708918486186e-6 1.8611855585014145e-6 2.3051551784149337e-6 2.7796930899776574e-6 3.263470692162403e-6 3.7303430704848267e-6 4.151490236519409e-6 4.498269337078975e-6 4.745402761905639e-6 4.8740153236803395e-6 4.8740153236803395e-6 4.745402761905639e-6 4.498269337078975e-6 4.151490236519405e-6 3.7303430704848267e-6 3.263470692162403e-6 2.7796930899776574e-6 2.3051551784149337e-6 1.8611855585014145e-6 1.4630708918486186e-6 1.1197659425728733e-6 8.344019823778828e-7 6.053542398020764e-7 4.275924597718886e-7 2.940604898402384e-7 1.9689264940429196e-7 1.2835373709777678e-7 8.146549888102353e-8 5.034137961561296e-8 3.0287448845822705e-8 1.7741341777273437e-8 1.0118040644152322e-8 5.618140103307805e-9 3.0372105473444726e-9; 1.6419398126916194e-9 3.0372105473444726e-9 5.4698920279304e-9 9.591108334634163e-9 1.6373632090899278e-8 2.7214944149646553e-8 4.404088682310447e-8 6.938903568369913e-8 1.0644170854928012e-7 1.5897140421510779e-7 2.3115983313045458e-7 3.272592438419348e-7 4.510842476340416e-7 6.053542398020764e-7 7.909475845251705e-7 1.0061714917927168e-6 1.2461849459797968e-6 1.502723858077153e-6 1.7642578193002676e-6 2.0166526840829715e-6 2.2443281409323587e-6 2.431799639052888e-6 2.5654019043369063e-6 2.6349308626683512e-6 2.6349308626683512e-6 2.5654019043369063e-6 2.431799639052888e-6 2.244328140932357e-6 2.0166526840829715e-6 1.7642578193002676e-6 1.502723858077153e-6 1.2461849459797968e-6 1.0061714917927168e-6 7.909475845251705e-7 6.053542398020754e-7 4.5108424763404237e-7 3.272592438419348e-7 2.3115983313045418e-7 1.5897140421510802e-7 1.0644170854928012e-7 6.938903568369913e-8 4.404088682310447e-8 2.72149441496466e-8 1.6373632090899278e-8 9.591108334634163e-9 5.4698920279304e-9 3.0372105473444726e-9 1.6419398126916194e-9])

To make the Gaussian Markov random field efficient we first precompute a bunch of quantities that allow us to scale things linearly with the number of image pixels. The returns a functional that accepts a single argument related to the correlation length of the field. The second argument defines the underlying random field of the Markov process. Here we are using a zero mean and unit variance Gaussian Markov random field. For this tutorial we will use the first order random field

julia
cprior = corr_image_prior(grid, dvis)
HierarchicalPrior(
	map: 
	ConditionalMarkov(
Random Field: VLBIImagePriors.GaussMarkovRandomField
Graph: MarkovRandomFieldGraph{1}(
dims: (48, 48)
)
)	hyper prior: 
	Truncated(Distributions.InverseGamma{Float64}(
invd: Distributions.Gamma{Float64}(α=1.0, θ=0.03623874772901475)
θ: 27.594772520225487
)
; lower=1.0, upper=96.0)

)

Putting everything together the total prior is then our image prior, a prior on the standard deviation of the MRF, and a prior on the fractional flux of the Gaussian component.

julia
prior = (
    c = cprior,
    σimg = truncated(Normal(0.0, 0.5); lower = 0.0),
    fg = Uniform(0.0, 1.0),
)
(c = HierarchicalPrior(
	map: 
	ConditionalMarkov(
Random Field: VLBIImagePriors.GaussMarkovRandomField
Graph: MarkovRandomFieldGraph{1}(
dims: (48, 48)
)
)	hyper prior: 
	Truncated(Distributions.InverseGamma{Float64}(
invd: Distributions.Gamma{Float64}(α=1.0, θ=0.03623874772901475)
θ: 27.594772520225487
)
; lower=1.0, upper=96.0)

)
, σimg = Truncated(Distributions.Normal{Float64}(μ=0.0, σ=0.5); lower=0.0), fg = Distributions.Uniform{Float64}(a=0.0, b=1.0))

Now we can construct our sky model.

julia
skym = SkyModel(sky, prior, grid; metadata = skymeta)
SkyModel
  with map: sky
   on grid: 
RectiGrid(
executor: ComradeBase.Serial()
Dimensions: 
(↓ X Sampled{Float64} LinRange{Float64}(-4.747133960864206e-10, 4.747133960864206e-10, 48) ForwardOrdered Regular Points,
→ Y Sampled{Float64} LinRange{Float64}(-4.747133960864206e-10, 4.747133960864206e-10, 48) ForwardOrdered Regular Points)
)
   )

Unlike other imaging examples (e.g., Imaging a Black Hole using only Closure Quantities) we also need to include a model for the instrument, i.e., gains. The gains will be broken into two components

  • Gain amplitudes which are typically known to 10-20%, except for LMT, which has amplitudes closer to 50-100%.

  • Gain phases which are more difficult to constrain and can shift rapidly.

julia
G = SingleStokesGain() do x
    lg = x.lg
    gp = x.gp
    return exp(lg + 1im * gp)
end

intpr = (
    lg = ArrayPrior(IIDSitePrior(ScanSeg(), Normal(0.0, 0.2)); LM = IIDSitePrior(ScanSeg(), Normal(0.0, 1.0))),
    gp = ArrayPrior(IIDSitePrior(ScanSeg(), DiagonalVonMises(0.0, inv(π^2))); refant = SEFDReference(0.0), phase = true),
)
intmodel = InstrumentModel(G, intpr)
InstrumentModel
  with Jones: SingleStokesGain
  with reference basis: PolarizedTypes.CirBasis()

To form the posterior we just combine the skymodel, instrument model and the data. To utilize gradients of the posterior we also need to load Enzyme.jl. Under the hood, Comrade will use Enzyme to compute the gradients of the posterior.

julia
using Enzyme
post = VLBIPosterior(skym, intmodel, dvis)
VLBIPosterior
ObservedSkyModel
  with map: sky
   on grid: 
FourierDualDomain(
Algorithm: VLBISkyModels.NFFTAlg{Float64, AbstractNFFTs.PrecomputeFlags, UInt32}(1, 1.0e-9, AbstractNFFTs.TENSOR, 0x00000000)
Image Domain: RectiGrid(
executor: ComradeBase.Serial()
Dimensions: 
(↓ X Sampled{Float64} LinRange{Float64}(-4.747133960864206e-10, 4.747133960864206e-10, 48) ForwardOrdered Regular Points,
→ Y Sampled{Float64} LinRange{Float64}(-4.747133960864206e-10, 4.747133960864206e-10, 48) ForwardOrdered Regular Points)
)
Visibility Domain: UnstructuredDomain(
executor: ComradeBase.Serial()
Dimensions: 
274-element StructArray(::Vector{Float64}, ::Vector{Float64}, ::Vector{Float64}, ::Vector{Float64}) with eltype @NamedTuple{U::Float64, V::Float64, Ti::Float64, Fr::Float64}:
 (U = -4.405690154666661e9, V = -4.523017159111106e9, Ti = 0.9166666567325592, Fr = 2.27070703125e11)
 (U = 787577.6145833326, V = -1.6838098888888871e6, Ti = 1.2166666388511658, Fr = 2.27070703125e11)
 (U = -4.444299918222218e9, V = -4.597825294222218e9, Ti = 1.2166666388511658, Fr = 2.27070703125e11)
 (U = 1.337045162666665e9, V = -3.765300401777774e9, Ti = 1.2166666388511658, Fr = 2.27070703125e11)
 (U = -1.336260540444443e9, V = 3.763616127999996e9, Ti = 1.2166666388511658, Fr = 2.27070703125e11)
 (U = 4.445088654222218e9, V = 4.596145080888884e9, Ti = 1.2166666388511658, Fr = 2.27070703125e11)
 (U = 5.781345607111105e9, V = 8.325259893333325e8, Ti = 1.2166666388511658, Fr = 2.27070703125e11)
 (U = 757554.6649305547, V = -1.6707483020833314e6, Ti = 1.516666665673256, Fr = 2.27070703125e11)
 (U = 1.4806382151111097e9, V = -3.741479615999996e9, Ti = 1.516666665673256, Fr = 2.27070703125e11)
 (U = -4.455366328888884e9, V = -4.673060451555551e9, Ti = 1.516666665673256, Fr = 2.27070703125e11)
 (U = -1.4798758791111097e9, V = 3.739809735111107e9, Ti = 1.516666665673256, Fr = 2.27070703125e11)
 (U = 4.456123861333328e9, V = 4.671391715555551e9, Ti = 1.516666665673256, Fr = 2.27070703125e11)
 (U = 5.936013027555549e9, V = 9.315912497777768e8, Ti = 1.516666665673256, Fr = 2.27070703125e11)
 (U = 722830.7065972214, V = -1.6582321493055536e6, Ti = 1.816666603088379, Fr = 2.27070703125e11)
 (U = -4.438811278222218e9, V = -4.748261176888884e9, Ti = 1.816666603088379, Fr = 2.27070703125e11)
 (U = 1.615060401777776e9, V = -3.715306943999996e9, Ti = 1.816666603088379, Fr = 2.27070703125e11)
 (U = -1.6143345528888872e9, V = 3.713649343999996e9, Ti = 1.816666603088379, Fr = 2.27070703125e11)
 (U = 4.439536184888884e9, V = 4.746598001777773e9, Ti = 1.816666603088379, Fr = 2.27070703125e11)
 (U = 6.053865329777771e9, V = 1.0329465084444433e9, Ti = 1.816666603088379, Fr = 2.27070703125e11)
 (U = 683620.5937499993, V = -1.6463396631944429e6, Ti = 2.1166666746139526, Fr = 2.27070703125e11)
 (U = -4.394744988444439e9, V = -4.822934485333328e9, Ti = 2.1166666746139526, Fr = 2.27070703125e11)
 (U = 1.7394626631111093e9, V = -3.686947740444441e9, Ti = 2.1166666746139526, Fr = 2.27070703125e11)
 (U = -1.7387777884444425e9, V = 3.6853017173333297e9, Ti = 2.1166666746139526, Fr = 2.27070703125e11)
 (U = 4.395427669333329e9, V = 4.821289443555551e9, Ti = 2.1166666746139526, Fr = 2.27070703125e11)
 (U = 6.134203007999993e9, V = 1.1359907128888876e9, Ti = 2.1166666746139526, Fr = 2.27070703125e11)
 (U = -1.8643712995555537e9, V = 3.651452543999996e9, Ti = 2.449999988079071, Fr = 2.27070703125e11)
 (U = 6.178850375111104e9, V = 1.2516714417777765e9, Ti = 2.449999988079071, Fr = 2.27070703125e11)
 (U = 4.314477112888885e9, V = 4.903116785777773e9, Ti = 2.449999988079071, Fr = 2.27070703125e11)
 (U = 587276.6180555549, V = -1.6236182569444429e6, Ti = 2.7500000596046448, Fr = 2.27070703125e11)
 (U = 1.9658264106666646e9, V = -3.6206990791111073e9, Ti = 2.7500000596046448, Fr = 2.27070703125e11)
 (U = -4.212778346666662e9, V = -4.976839239111106e9, Ti = 2.7500000596046448, Fr = 2.27070703125e11)
 (U = -1.965238759111109e9, V = 3.6190756053333297e9, Ti = 2.7500000596046448, Fr = 2.27070703125e11)
 (U = 6.178607487999993e9, V = 1.356139594666665e9, Ti = 2.7500000596046448, Fr = 2.27070703125e11)
 (U = 4.213364039111107e9, V = 4.975216568888884e9, Ti = 2.7500000596046448, Fr = 2.27070703125e11)
 (U = 535806.2881944438, V = -1.6141240694444429e6, Ti = 3.0500001311302185, Fr = 2.27070703125e11)
 (U = 2.0544588266666646e9, V = -3.586710200888885e9, Ti = 3.0500001311302185, Fr = 2.27070703125e11)
 (U = -4.085601863111107e9, V = -5.046995640888884e9, Ti = 3.0500001311302185, Fr = 2.27070703125e11)
 (U = -2.0539194026666646e9, V = 3.5850976213333297e9, Ti = 3.0500001311302185, Fr = 2.27070703125e11)
 (U = 4.0861423288888845e9, V = 5.045379299555551e9, Ti = 3.0500001311302185, Fr = 2.27070703125e11)
 (U = 6.140063729777771e9, V = 1.4602817422222207e9, Ti = 3.0500001311302185, Fr = 2.27070703125e11)
 (U = 481010.2230902773, V = -1.6055275520833316e6, Ti = 3.3499998450279236, Fr = 2.27070703125e11)
 (U = 2.130349845333331e9, V = -3.5513328995555515e9, Ti = 3.3499998450279236, Fr = 2.27070703125e11)
 (U = -3.933104519111107e9, V = -5.114784767999995e9, Ti = 3.3499998450279236, Fr = 2.27070703125e11)
 (U = -2.1298682808888867e9, V = 3.5497276302222185e9, Ti = 3.3499998450279236, Fr = 2.27070703125e11)
 (U = 6.063457521777771e9, V = 1.5634511359999983e9, Ti = 3.3499998450279236, Fr = 2.27070703125e11)
 (U = 3.9335861404444404e9, V = 5.113178993777773e9, Ti = 3.3499998450279236, Fr = 2.27070703125e11)
 (U = 416648.93055555515, V = -1.5970948472222206e6, Ti = 3.6833333373069763, Fr = 2.27070703125e11)
 (U = 2.578938929777775e9, V = -4.734788238222218e9, Ti = 3.6833333373069763, Fr = 2.27070703125e11)
 (U = -3.735105663999996e9, V = -5.186825713777772e9, Ti = 3.6833333373069763, Fr = 2.27070703125e11)
 (U = 2.1991734328888865e9, V = -3.5106569031111073e9, Ti = 3.6833333373069763, Fr = 2.27070703125e11)
 (U = -2.5785191964444413e9, V = 4.733192191999995e9, Ti = 3.6833333373069763, Fr = 2.27070703125e11)
 (U = -2.198756124444442e9, V = 3.509060273777774e9, Ti = 3.6833333373069763, Fr = 2.27070703125e11)
 (U = 3.797641537777774e8, V = -1.2241308906666653e9, Ti = 3.6833333373069763, Fr = 2.27070703125e11)
 (U = 3.735526620444441e9, V = 5.185227192888883e9, Ti = 3.6833333373069763, Fr = 2.27070703125e11)
 (U = 6.314048611555549e9, V = 4.520333351111106e8, Ti = 3.6833333373069763, Fr = 2.27070703125e11)
 (U = 5.934280931555549e9, V = 1.6761668835555537e9, Ti = 3.6833333373069763, Fr = 2.27070703125e11)
 (U = 364284.04464285675, V = -1.591376535714284e6, Ti = 3.98333340883255, Fr = 2.27070703125e11)
 (U = -3.5324296035555515e9, V = -5.248264732444439e9, Ti = 3.98333340883255, Fr = 2.27070703125e11)
 (U = 2.246733795555553e9, V = -3.4730721351111073e9, Ti = 3.98333340883255, Fr = 2.27070703125e11)
 (U = 2.7040352142222195e9, V = -4.690125752888884e9, Ti = 3.98333340883255, Fr = 2.27070703125e11)
 (U = -2.687901513142854e9, V = 4.694663167999995e9, Ti = 3.98333340883255, Fr = 2.27070703125e11)
 (U = -2.2408562468571405e9, V = 3.476569380571425e9, Ti = 3.98333340883255, Fr = 2.27070703125e11)
 (U = 4.573011991111106e8, V = -1.2170540408888876e9, Ti = 3.98333340883255, Fr = 2.27070703125e11)
 (U = 3.561408036571425e9, V = 5.238641298285708e9, Ti = 3.98333340883255, Fr = 2.27070703125e11)
 (U = 5.779150734222217e9, V = 1.7751998862222204e9, Ti = 3.98333340883255, Fr = 2.27070703125e11)
 (U = 6.236464440888882e9, V = 5.581377528888882e8, Ti = 3.98333340883255, Fr = 2.27070703125e11)
 (U = 292810.0842803027, V = -1.5850529772727257e6, Ti = 4.283333480358124, Fr = 2.27070703125e11)
 (U = 2.81236712533333e9, V = -4.643490716444439e9, Ti = 4.283333480358124, Fr = 2.27070703125e11)
 (U = 2.280368668444442e9, V = -3.43479938133333e9, Ti = 4.283333480358124, Fr = 2.27070703125e11)
 (U = -3.307854179555552e9, V = -5.306092316444439e9, Ti = 4.283333480358124, Fr = 2.27070703125e11)
 (U = -2.812399111529409e9, V = 4.641748931764702e9, Ti = 4.283333480358124, Fr = 2.27070703125e11)
 (U = -2.280159088941174e9, V = 3.43308967905882e9, Ti = 4.283333480358124, Fr = 2.27070703125e11)
 (U = 5.3200019644444394e8, V = -1.2086904924444432e9, Ti = 4.283333480358124, Fr = 2.27070703125e11)
 (U = 3.307383679999996e9, V = 5.304687194352936e9, Ti = 4.283333480358124, Fr = 2.27070703125e11)
 (U = 6.120220771555549e9, V = 6.626024835555549e8, Ti = 4.283333480358124, Fr = 2.27070703125e11)
 (U = 5.588211185777772e9, V = 1.8712975182222202e9, Ti = 4.283333480358124, Fr = 2.27070703125e11)
 (U = 238138.85156249977, V = -1.5812476874999984e6, Ti = 4.583333194255829, Fr = 2.27070703125e11)
 (U = 5.092881095111105e9, V = -4.199598421333329e9, Ti = 4.583333194255829, Fr = 2.27070703125e11)
 (U = 2.903268209777775e9, V = -4.595170062222218e9, Ti = 4.583333194255829, Fr = 2.27070703125e11)
 (U = 2.299867178666664e9, V = -3.396077226666663e9, Ti = 4.583333194255829, Fr = 2.27070703125e11)
 (U = -3.062767004444441e9, V = -5.359950862222217e9, Ti = 4.583333194255829, Fr = 2.27070703125e11)
 (U = 5.050181503999994e9, V = -4.210810559999995e9, Ti = 4.583333194255829, Fr = 2.27070703125e11)
 (U = -2.8907034879999967e9, V = 4.600886271999995e9, Ti = 4.583333194255829, Fr = 2.27070703125e11)
 (U = 2.189606492444442e9, V = 3.955690115555551e8, Ti = 4.583333194255829, Fr = 2.27070703125e11)
 (U = -2.2976930559999976e9, V = 3.400280575999996e9, Ti = 4.583333194255829, Fr = 2.27070703125e11)
 (U = 2.79301371733333e9, V = -8.035218719999992e8, Ti = 4.583333194255829, Fr = 2.27070703125e11)
 (U = 6.034046399999994e8, V = -1.199091349333332e9, Ti = 4.583333194255829, Fr = 2.27070703125e11)
 (U = 3.1008046079999967e9, V = 5.350639359999994e9, Ti = 4.583333194255829, Fr = 2.27070703125e11)
 (U = 8.155639423999991e9, V = 1.1603480675555544e9, Ti = 4.583333194255829, Fr = 2.27070703125e11)
 (U = 5.966032455111105e9, V = 7.64783736888888e8, Ti = 4.583333194255829, Fr = 2.27070703125e11)
 (U = 5.36262916266666e9, V = 1.963875015111109e9, Ti = 4.583333194255829, Fr = 2.27070703125e11)
 (U = 163232.87499999983, V = -1.5773965972222206e6, Ti = 4.916666686534882, Fr = 2.27070703125e11)
 (U = 5.388079544888883e9, V = -4.101135303111107e9, Ti = 4.916666686534882, Fr = 2.27070703125e11)
 (U = 2.3048084479999976e9, V = -3.3528182257777743e9, Ti = 4.916666686534882, Fr = 2.27070703125e11)
 (U = -2.768253631999997e9, V = -5.414730879999994e9, Ti = 4.916666686534882, Fr = 2.27070703125e11)
 (U = 5.388055679999994e9, V = -4.101249678222218e9, Ti = 4.916666686534882, Fr = 2.27070703125e11)
 (U = 2.9831417315555525e9, V = -4.539868501333328e9, Ti = 4.916666686534882, Fr = 2.27070703125e11)
 (U = 5.356877311999994e9, V = -4.110922695111107e9, Ti = 4.916666686534882, Fr = 2.27070703125e11)
 (U = -2.975202047999997e9, V = 4.544589937777773e9, Ti = 4.916666686534882, Fr = 2.27070703125e11)
 (U = 2.4049356302222195e9, V = 4.387324826666662e8, Ti = 4.916666686534882, Fr = 2.27070703125e11)
 (U = 2.4049011413333306e9, V = 4.3861474933333284e8, Ti = 4.916666686534882, Fr = 2.27070703125e11)
 (U = -2.305032220444442e9, V = 3.356104334222219e9, Ti = 4.916666686534882, Fr = 2.27070703125e11)
 (U = 3.083277411555552e9, V = -7.48315953777777e8, Ti = 4.916666686534882, Fr = 2.27070703125e11)
 (U = 3.0832264319999967e9, V = -7.484370453333325e8, Ti = 4.916666686534882, Fr = 2.27070703125e11)
 (U = 6.783317208888881e8, V = -1.187050453333332e9, Ti = 4.916666686534882, Fr = 2.27070703125e11)
 (U = 2.8027259164444413e9, V = 5.407289457777772e9, Ti = 4.916666686534882, Fr = 2.27070703125e11)
 (U = 8.156328149333324e9, V = 1.3135985066666653e9, Ti = 4.916666686534882, Fr = 2.27070703125e11)
 (U = 8.156297400888881e9, V = 1.3134812622222207e9, Ti = 4.916666686534882, Fr = 2.27070703125e11)
 (U = 5.751382712888883e9, V = 8.74867279999999e8, Ti = 4.916666686534882, Fr = 2.27070703125e11)
 (U = 5.073054606222218e9, V = 2.0619144959999979e9, Ti = 4.916666686534882, Fr = 2.27070703125e11)
 (U = -5.356795278222217e9, V = 4.1110577777777734e9, Ti = 4.916666686534882, Fr = 2.27070703125e11)
 (U = 30750.899088541635, V = 116905.95138888876, Ti = 4.916666686534882, Fr = 2.27070703125e11)
 (U = 99881.10751488084, V = -1.575297547619046e6, Ti = 5.183333337306976, Fr = 2.27070703125e11)
 (U = 5.59465437866666e9, V = -4.018611114666662e9, Ti = 5.183333337306976, Fr = 2.27070703125e11)
 (U = 3.030624810666663e9, V = -4.494683207111106e9, Ti = 5.183333337306976, Fr = 2.27070703125e11)
 (U = 2.2960602239999976e9, V = -3.3182473599999967e9, Ti = 5.183333337306976, Fr = 2.27070703125e11)
 (U = 5.594626830222217e9, V = -4.0187286115555515e9, Ti = 5.183333337306976, Fr = 2.27070703125e11)
 (U = -2.5172720142222195e9, V = -5.45444669155555e9, Ti = 5.183333337306976, Fr = 2.27070703125e11)
 (U = 5.579493814857137e9, V = -4.023605065142853e9, Ti = 5.183333337306976, Fr = 2.27070703125e11)
 (U = -3.0274068479999967e9, V = 4.496662674285709e9, Ti = 5.183333337306976, Fr = 2.27070703125e11)
 (U = 2.5640353279999976e9, V = 4.7607165511111057e8, Ti = 5.183333337306976, Fr = 2.27070703125e11)
 (U = 2.564000689777775e9, V = 4.7595243111111057e8, Ti = 5.183333337306976, Fr = 2.27070703125e11)
 (U = -2.2970926933333306e9, V = 3.3193651078095202e9, Ti = 5.183333337306976, Fr = 2.27070703125e11)
 (U = 3.2985927679999967e9, V = -7.003643235555549e8, Ti = 5.183333337306976, Fr = 2.27070703125e11)
 (U = 7.345665155555549e8, V = -1.1764336924444432e9, Ti = 5.183333337306976, Fr = 2.27070703125e11)
 (U = 3.2985603128888855e9, V = -7.004829137777771e8, Ti = 5.183333337306976, Fr = 2.27070703125e11)
 (U = 2.5374800944761877e9, V = 5.44992148723809e9, Ti = 5.183333337306976, Fr = 2.27070703125e11)
 (U = 8.111924807111103e9, V = 1.4358354488888874e9, Ti = 5.183333337306976, Fr = 2.27070703125e11)
 (U = 4.813340984888884e9, V = 2.136196949333331e9, Ti = 5.183333337306976, Fr = 2.27070703125e11)
 (U = 5.547893447111105e9, V = 9.597677279999989e8, Ti = 5.183333337306976, Fr = 2.27070703125e11)
 (U = 8.111898709333324e9, V = 1.4357162488888874e9, Ti = 5.183333337306976, Fr = 2.27070703125e11)
 (U = -5.579455049142851e9, V = 4.0237272990476146e9, Ti = 5.183333337306976, Fr = 2.27070703125e11)
 (U = 27098.63357204858, V = 117340.62348090265, Ti = 5.183333337306976, Fr = 2.27070703125e11)
 (U = 45690.53906249995, V = -1.5743099999999984e6, Ti = 5.449999988079071, Fr = 2.27070703125e11)
 (U = 5.773825464888883e9, V = -3.933187214222218e9, Ti = 5.449999988079071, Fr = 2.27070703125e11)
 (U = 2.2760630257777753e9, V = -3.283890659555552e9, Ti = 5.449999988079071, Fr = 2.27070703125e11)
 (U = 3.0632664319999967e9, V = -4.448890581333328e9, Ti = 5.449999988079071, Fr = 2.27070703125e11)
 (U = 5.773791857777772e9, V = -3.9333102648888845e9, Ti = 5.449999988079071, Fr = 2.27070703125e11)
 (U = -2.253957283555553e9, V = -5.490298367999994e9, Ti = 5.449999988079071, Fr = 2.27070703125e11)
 (U = 5.744240127999994e9, V = -3.946954495999996e9, Ti = 5.449999988079071, Fr = 2.27070703125e11)
 (U = -3.0586252799999967e9, V = 4.455449599999995e9, Ti = 5.449999988079071, Fr = 2.27070703125e11)
 (U = 2.710558015999997e9, V = 5.1570311111111057e8, Ti = 5.449999988079071, Fr = 2.27070703125e11)
 (U = 2.7105363982222195e9, V = 5.1558591999999946e8, Ti = 5.449999988079071, Fr = 2.27070703125e11)
 (U = -2.2804378879999976e9, V = 3.2883723519999967e9, Ti = 5.449999988079071, Fr = 2.27070703125e11)
 (U = 3.4977656675555515e9, V = -6.492943324444438e8, Ti = 5.449999988079071, Fr = 2.27070703125e11)
 (U = 7.872002471111102e8, V = -1.164999438222221e9, Ti = 5.449999988079071, Fr = 2.27070703125e11)
 (U = 3.497740145777774e9, V = -6.494127164444437e8, Ti = 5.449999988079071, Fr = 2.27070703125e11)
 (U = 2.3015257599999976e9, V = 5.482690303999994e9, Ti = 5.449999988079071, Fr = 2.27070703125e11)
 (U = 8.027784391111103e9, V = 1.5571086684444427e9, Ti = 5.449999988079071, Fr = 2.27070703125e11)
 (U = 5.317234005333327e9, V = 1.0414041991111101e9, Ti = 5.449999988079071, Fr = 2.27070703125e11)
 (U = 8.027758634666658e9, V = 1.556993443555554e9, Ti = 5.449999988079071, Fr = 2.27070703125e11)
 (U = 4.530021788444439e9, V = 2.206406378666664e9, Ti = 5.449999988079071, Fr = 2.27070703125e11)
 (U = -5.744238079999993e9, V = 3.9470612479999957e9, Ti = 5.449999988079071, Fr = 2.27070703125e11)
 (U = 23313.64084201386, V = 117719.41384548598, Ti = 5.449999988079071, Fr = 2.27070703125e11)
 (U = 5.924709319111105e9, V = -3.8452848639999957e9, Ti = 5.716666638851166, Fr = 2.27070703125e11)
 (U = 5.924688497777772e9, V = -3.845403854222218e9, Ti = 5.716666638851166, Fr = 2.27070703125e11)
 (U = 3.08089890133333e9, V = -4.402725802666661e9, Ti = 5.716666638851166, Fr = 2.27070703125e11)
 (U = 2.244918698666664e9, V = -3.2499217208888855e9, Ti = 5.716666638851166, Fr = 2.27070703125e11)
 (U = 2.8438056746666636e9, V = 5.57437258666666e8, Ti = 5.716666638851166, Fr = 2.27070703125e11)
 (U = 2.8437909475555525e9, V = 5.573208195555549e8, Ti = 5.716666638851166, Fr = 2.27070703125e11)
 (U = 3.6797841706666627e9, V = -5.953657653333327e8, Ti = 5.716666638851166, Fr = 2.27070703125e11)
 (U = 3.679770524444441e9, V = -5.954820355555549e8, Ti = 5.716666638851166, Fr = 2.27070703125e11)
 (U = 8.35980113777777e8, V = -1.1528026026666653e9, Ti = 5.716666638851166, Fr = 2.27070703125e11)
 (U = 19414.168077256923, V = 118040.4856770832, Ti = 5.716666638851166, Fr = 2.27070703125e11)
 (U = -92093.64322916657, V = -1.5751088611111094e6, Ti = 6.049999952316284, Fr = 2.27070703125e11)
 (U = 3.081715143111108e9, V = -4.344833365333329e9, Ti = 6.049999952316284, Fr = 2.27070703125e11)
 (U = 6.07242564266666e9, V = -3.732578396444441e9, Ti = 6.049999952316284, Fr = 2.27070703125e11)
 (U = 2.190544647111109e9, V = -3.208251676444441e9, Ti = 6.049999952316284, Fr = 2.27070703125e11)
 (U = 6.05836305066666e9, V = -3.7438256071111073e9, Ti = 6.049999952316284, Fr = 2.27070703125e11)
 (U = -3.0829770239999967e9, V = 4.34976796444444e9, Ti = 6.049999952316284, Fr = 2.27070703125e11)
 (U = 2.99071138133333e9, V = 6.122509119999994e8, Ti = 6.049999952316284, Fr = 2.27070703125e11)
 (U = -2.1976609848888865e9, V = 3.211308942222219e9, Ti = 6.049999952316284, Fr = 2.27070703125e11)
 (U = 3.881882190222218e9, V = -5.2432489066666615e8, Ti = 6.049999952316284, Fr = 2.27070703125e11)
 (U = 8.911678026666657e8, V = -1.1365769386666653e9, Ti = 6.049999952316284, Fr = 2.27070703125e11)
 (U = -150875.45052083317, V = -1.5769051249999984e6, Ti = 6.316666603088379, Fr = 2.27070703125e11)
 (U = 6.157184056888882e9, V = -3.6406887111111073e9, Ti = 6.316666603088379, Fr = 2.27070703125e11)
 (U = 6.157173617777771e9, V = -3.6408075306666627e9, Ti = 6.316666603088379, Fr = 2.27070703125e11)
 (U = 2.1349524657777758e9, V = -3.175750087111108e9, Ti = 6.316666603088379, Fr = 2.27070703125e11)
 (U = 3.065380067555552e9, V = -4.298637105777774e9, Ti = 6.316666603088379, Fr = 2.27070703125e11)
 (U = 6.14621098666666e9, V = -3.6535959466666627e9, Ti = 6.316666603088379, Fr = 2.27070703125e11)
 (U = -3.0691708159999967e9, V = 4.304277162666662e9, Ti = 6.316666603088379, Fr = 2.27070703125e11)
 (U = 3.091809799111108e9, V = 6.579570684444437e8, Ti = 6.316666603088379, Fr = 2.27070703125e11)
 (U = 3.091795150222219e9, V = 6.578363697777771e8, Ti = 6.316666603088379, Fr = 2.27070703125e11)
 (U = -2.1445705386666644e9, V = 3.1792032426666636e9, Ti = 6.316666603088379, Fr = 2.27070703125e11)
 (U = 4.0222341617777734e9, V = -4.6493393599999946e8, Ti = 6.316666603088379, Fr = 2.27070703125e11)
 (U = 4.0222294044444404e9, V = -4.650499119999995e8, Ti = 6.316666603088379, Fr = 2.27070703125e11)
 (U = 9.304277457777768e8, V = -1.1228893119999988e9, Ti = 6.316666603088379, Fr = 2.27070703125e11)
 (U = -6.146197247999993e9, V = 3.653717759999996e9, Ti = 6.316666603088379, Fr = 2.27070703125e11)
 (U = 10322.971652560753, V = 118544.06944444432, Ti = 6.316666603088379, Fr = 2.27070703125e11)
 (U = 6.211784632888882e9, V = -3.5477483875555515e9, Ti = 6.583333253860474, Fr = 2.27070703125e11)
 (U = 6.211780195555549e9, V = -3.547863480888885e9, Ti = 6.583333253860474, Fr = 2.27070703125e11)
 (U = 2.068902503111109e9, V = -3.1441626808888855e9, Ti = 6.583333253860474, Fr = 2.27070703125e11)
 (U = 3.034031331555552e9, V = -4.2528075377777734e9, Ti = 6.583333253860474, Fr = 2.27070703125e11)
 (U = 3.177752227555552e9, V = 7.050634364444437e8, Ti = 6.583333253860474, Fr = 2.27070703125e11)
 (U = 3.177752810666663e9, V = 7.049488924444437e8, Ti = 6.583333253860474, Fr = 2.27070703125e11)
 (U = 4.1428866915555515e9, V = -4.035797448888885e8, Ti = 6.583333253860474, Fr = 2.27070703125e11)
 (U = 9.651305439999989e8, V = -1.108645550222221e9, Ti = 6.583333253860474, Fr = 2.27070703125e11)
 (U = 4.1428811519999957e9, V = -4.036979662222218e8, Ti = 6.583333253860474, Fr = 2.27070703125e11)
 (U = 6186.6973334418335, V = 118668.12304687487, Ti = 6.583333253860474, Fr = 2.27070703125e11)
 (U = -266936.3645833331, V = -1.583141166666665e6, Ti = 6.849999904632568, Fr = 2.27070703125e11)
 (U = 6.235954645333326e9, V = -3.4543274026666627e9, Ti = 6.849999904632568, Fr = 2.27070703125e11)
 (U = 2.987818652444441e9, V = -4.207558826666662e9, Ti = 6.849999904632568, Fr = 2.27070703125e11)
 (U = 6.235956423111104e9, V = -3.454211199999996e9, Ti = 6.849999904632568, Fr = 2.27070703125e11)
 (U = 1.9927185599999979e9, V = -3.1136442239999967e9, Ti = 6.849999904632568, Fr = 2.27070703125e11)
 (U = 6.23440708266666e9, V = -3.468812543999996e9, Ti = 6.849999904632568, Fr = 2.27070703125e11)
 (U = -2.9972223146666636e9, V = 4.213749674666662e9, Ti = 6.849999904632568, Fr = 2.27070703125e11)
 (U = 3.2481352177777743e9, V = 7.533482702222215e8, Ti = 6.849999904632568, Fr = 2.27070703125e11)
 (U = 3.2481338097777743e9, V = 7.532299875555547e8, Ti = 6.849999904632568, Fr = 2.27070703125e11)
 (U = -2.006915669333331e9, V = 3.1172445866666636e9, Ti = 6.849999904632568, Fr = 2.27070703125e11)
 (U = 4.243231630222218e9, V = -3.4069073333333296e8, Ti = 6.849999904632568, Fr = 2.27070703125e11)
 (U = 4.243239466666662e9, V = -3.4056780622222185e8, Ti = 6.849999904632568, Fr = 2.27070703125e11)
 (U = 9.951025973333323e8, V = -1.0939171839999988e9, Ti = 6.849999904632568, Fr = 2.27070703125e11)
 (U = -6.23440452266666e9, V = 3.4689294506666627e9, Ti = 6.849999904632568, Fr = 2.27070703125e11)
 (U = 2020.1427137586786, V = 118729.78884548598, Ti = 6.849999904632568, Fr = 2.27070703125e11)
 (U = -343950.03624999966, V = -1.589413219999998e6, Ti = 7.183333396911621, Fr = 2.27070703125e11)
 (U = 6.223212003555549e9, V = -3.337163527111108e9, Ti = 7.183333396911621, Fr = 2.27070703125e11)
 (U = 6.223215416888882e9, V = -3.337282766222219e9, Ti = 7.183333396911621, Fr = 2.27070703125e11)
 (U = 2.909506687999997e9, V = -4.152156949333329e9, Ti = 7.183333396911621, Fr = 2.27070703125e11)
 (U = 1.8838133048888867e9, V = -3.07722722133333e9, Ti = 7.183333396911621, Fr = 2.27070703125e11)
 (U = 6.225276682239994e9, V = -3.340931624959996e9, Ti = 7.183333396911621, Fr = 2.27070703125e11)
 (U = -2.913978071039997e9, V = 4.1530716876799955e9, Ti = 7.183333396911621, Fr = 2.27070703125e11)
 (U = 3.313704135111108e9, V = 8.14994561777777e8, Ti = 7.183333396911621, Fr = 2.27070703125e11)
 (U = 3.3137085368888855e9, V = 8.148773351111102e8, Ti = 7.183333396911621, Fr = 2.27070703125e11)
 (U = -1.889491799039998e9, V = 3.077254901759997e9, Ti = 7.183333396911621, Fr = 2.27070703125e11)
 (U = 4.339403889777773e9, V = -2.5993809599999973e8, Ti = 7.183333396911621, Fr = 2.27070703125e11)
 (U = 4.339408696888885e9, V = -2.6005522088888863e8, Ti = 7.183333396911621, Fr = 2.27070703125e11)
 (U = 1.0256997831111101e9, V = -1.0749325617777767e9, Ti = 7.183333396911621, Fr = 2.27070703125e11)
 (U = -6.225278402559994e9, V = 3.3410462617599964e9, Ti = 7.183333396911621, Fr = 2.27070703125e11)
 (U = -3199.665127224389, V = 118718.70724826376, Ti = 7.183333396911621, Fr = 2.27070703125e11)
 (U = -401459.5972222218, V = -1.5953378194444429e6, Ti = 7.450000047683716, Fr = 2.27070703125e11)
 (U = 6.178709461333326e9, V = -3.243972359111108e9, Ti = 7.450000047683716, Fr = 2.27070703125e11)
 (U = 2.8307861475555525e9, V = -4.1090220657777734e9, Ti = 7.450000047683716, Fr = 2.27070703125e11)
 (U = 6.178720483555549e9, V = -3.244096682666663e9, Ti = 7.450000047683716, Fr = 2.27070703125e11)
 (U = 1.786251743999998e9, V = -3.049647466666663e9, Ti = 7.450000047683716, Fr = 2.27070703125e11)
 (U = 6.179113016888882e9, V = -3.242380138666663e9, Ti = 7.450000047683716, Fr = 2.27070703125e11)
 (U = 3.347918890666663e9, V = 8.650494186666657e8, Ti = 7.450000047683716, Fr = 2.27070703125e11)
 (U = -2.8311915662222195e9, V = 4.107428693333329e9, Ti = 7.450000047683716, Fr = 2.27070703125e11)
 (U = 3.347926620444441e9, V = 8.649315555555546e8, Ti = 7.450000047683716, Fr = 2.27070703125e11)
 (U = -1.7866596017777758e9, V = 3.048053795555552e9, Ti = 7.450000047683716, Fr = 2.27070703125e11)
 (U = 4.392453788444439e9, V = -1.9432878977777755e8, Ti = 7.450000047683716, Fr = 2.27070703125e11)
 (U = 4.392462890666661e9, V = -1.944448759999998e8, Ti = 7.450000047683716, Fr = 2.27070703125e11)
 (U = 1.0445360586666656e9, V = -1.0593764586666656e9, Ti = 7.450000047683716, Fr = 2.27070703125e11)
 (U = -6.179119544888882e9, V = 3.2424975217777743e9, Ti = 7.450000047683716, Fr = 2.27070703125e11)
 (U = -7359.556450737839, V = 118639.36892361098, Ti = 7.450000047683716, Fr = 2.27070703125e11)
 (U = -453840.1640624995, V = -1.6017645034722206e6, Ti = 7.7166666984558105, Fr = 2.27070703125e11)
 (U = 2.7382116835555525e9, V = -4.0671804444444404e9, Ti = 7.7166666984558105, Fr = 2.27070703125e11)
 (U = 6.103944163555549e9, V = -3.151683946666663e9, Ti = 7.7166666984558105, Fr = 2.27070703125e11)
 (U = 1.6799544319999983e9, V = -3.023603527111108e9, Ti = 7.7166666984558105, Fr = 2.27070703125e11)
 (U = 6.103956750222216e9, V = -3.1518036479999967e9, Ti = 7.7166666984558105, Fr = 2.27070703125e11)
 (U = 6.10439825066666e9, V = -3.150082410666663e9, Ti = 7.7166666984558105, Fr = 2.27070703125e11)
 (U = 3.3657344497777743e9, V = 9.154982524444435e8, Ti = 7.7166666984558105, Fr = 2.27070703125e11)
 (U = -2.7386646755555525e9, V = 4.065578346666662e9, Ti = 7.7166666984558105, Fr = 2.27070703125e11)
 (U = 3.365744974222219e9, V = 9.153749564444435e8, Ti = 7.7166666984558105, Fr = 2.27070703125e11)
 (U = -1.6804040142222204e9, V = 3.02200078933333e9, Ti = 7.7166666984558105, Fr = 2.27070703125e11)
 (U = 4.423988807111106e9, V = -1.2808080399999987e8, Ti = 7.7166666984558105, Fr = 2.27070703125e11)
 (U = 4.424001663999995e9, V = -1.2819521311111099e8, Ti = 7.7166666984558105, Fr = 2.27070703125e11)
 (U = 1.0582553617777767e9, V = -1.0435755999999989e9, Ti = 7.7166666984558105, Fr = 2.27070703125e11)
 (U = -6.104410581333326e9, V = 3.150201891555552e9, Ti = 7.7166666984558105, Fr = 2.27070703125e11)
 (U = -11483.574490017349, V = 118497.78168402765, Ti = 7.7166666984558105, Fr = 2.27070703125e11)
 (U = -504000.3645833328, V = -1.6089620833333316e6, Ti = 7.983333349227905, Fr = 2.27070703125e11)
 (U = 5.999277468444439e9, V = -3.060741980444441e9, Ti = 7.983333349227905, Fr = 2.27070703125e11)
 (U = 2.6322108302222195e9, V = -4.0268253795555515e9, Ti = 7.983333349227905, Fr = 2.27070703125e11)
 (U = 5.999288348444438e9, V = -3.060856782222219e9, Ti = 7.983333349227905, Fr = 2.27070703125e11)
 (U = 1.565423818666665e9, V = -2.999217870222219e9, Ti = 7.983333349227905, Fr = 2.27070703125e11)
 (U = 5.999780252444438e9, V = -3.0591319537777743e9, Ti = 7.983333349227905, Fr = 2.27070703125e11)
 (U = -2.6327225742222195e9, V = 4.025219150222218e9, Ti = 7.983333349227905, Fr = 2.27070703125e11)
 (U = 3.3670615039999967e9, V = 9.660856764444433e8, Ti = 7.983333349227905, Fr = 2.27070703125e11)
 (U = 3.367077162666663e9, V = 9.659668408888879e8, Ti = 7.983333349227905, Fr = 2.27070703125e11)
 (U = 4.433867007999995e9, V = -6.164012344444438e7, Ti = 7.983333349227905, Fr = 2.27070703125e11)
 (U = -1.5659301759999983e9, V = 2.99760936533333e9, Ti = 7.983333349227905, Fr = 2.27070703125e11)
 (U = 4.433851406222218e9, V = -6.1523151111111045e7, Ti = 7.983333349227905, Fr = 2.27070703125e11)
 (U = 1.0667900266666656e9, V = -1.0276084337777767e9, Ti = 7.983333349227905, Fr = 2.27070703125e11)
 (U = -5.999801315555549e9, V = 3.059254300444441e9, Ti = 7.983333349227905, Fr = 2.27070703125e11)
 (U = -15551.297851562484, V = 118294.64453124987, Ti = 7.983333349227905, Fr = 2.27070703125e11)
)
   )
ObservedInstrumentModel
  with Jones: SingleStokesGain
  with reference basis: PolarizedTypes.CirBasis()Data Products: Comrade.EHTVisibilityDatum

Optimization and Sampling

Now we need to actually compute our image. For this we will first follow standard approaches in VLBI and find the maximum a posteriori (MAP) estimate of the image and instrument model. This is done using the comrade_opt function which accepts the posterior, an optimization algorithm, and some keyword arguments. For this tutorial we will use the L-BFGS algorithm. For more information about the optimization algorithms available see the Optimization.jl docs.

julia
using Optimization, OptimizationLBFGSB
xopt, sol = comrade_opt(
    post, LBFGSB(); initial_params = prior_sample(rng, post),
    maxiters = 2000, g_tol = 1.0e-1
);

Warning

Fitting gains tends to be very difficult, meaning that optimization can take a lot longer. The upside is that we usually get nicer images.

First we will evaluate our fit by plotting the residuals

julia
using CairoMakie
using DisplayAs
res = residuals(post, xopt)
plotfields(res[1], :uvdist, :res) |> DisplayAs.PNG |> DisplayAs.Text

These look reasonable, although there may be some minor overfitting. This could be improved in a few ways, but that is beyond the goal of this quick tutorial. Plotting the image, we see that we have a much cleaner version of the closure-only image from Imaging a Black Hole using only Closure Quantities.

julia
g = imagepixels(fovx, fovy, 128, 128)
img = intensitymap(skymodel(post, xopt), g)
imageviz(img, size = (500, 400)) |> DisplayAs.PNG |> DisplayAs.Text

Because we also fit the instrument model, we can inspect their parameters. First, let's query the posterior object with the optimal parameters to get the instrument model.

julia
intopt = instrumentmodel(post, xopt)
126-element Comrade.SiteArray{ComplexF64, 1, Vector{ComplexF64}, Vector{Comrade.IntegrationTime{Int64, Float64}}, Vector{Comrade.FrequencyChannel{Float64, Int64}}, Vector{Symbol}}:
    0.9998015073462388 + 0.0im
   -0.5581776916218577 + 0.830193972522752im
    1.0086638657310052 + 0.0im
   -0.5979527964582398 - 0.8356355568315927im
    -0.384944313544361 + 0.7307230296510537im
   -0.4497493705150899 + 0.9878814238966094im
     1.004540506884221 + 0.0im
   -0.6400645579612088 - 0.8070653837364133im
   -0.5230943172903714 + 0.8271236375719124im
  -0.22803162373290542 + 1.0165686547852135im
    1.0191408193405431 + 0.0im
    -0.660257916075897 - 0.77685076579337im
   -0.5250564597742826 + 0.7494170070419572im
   -0.0396539810003907 + 1.014429554635044im
    0.9969748443822941 + 0.0im
   -0.7128007440090385 - 0.7537049355897725im
   -0.5822441182359516 + 0.7289911726790523im
   0.14427458213199337 + 0.9489247130843392im
    1.0399535778052285 + 0.0im
   -0.6277345895403884 + 0.2624734344372759im
   -0.7841502501797877 - 0.40420599330800666im
    1.0210793627532135 + 0.0im
   -0.7533683535269438 - 0.6868979916192778im
   -0.3942925108816618 + 0.4878195973630148im
    0.4948510435862709 + 0.8492730276475869im
     1.009343174617229 + 0.0im
   -0.7750813681938435 - 0.6621454399314219im
   -0.3671298016107702 + 0.4581538317245425im
    0.6661262956570515 + 0.75640768315302im
    1.0210080798251113 + 0.0im
   -0.7943862942089152 - 0.637667951654791im
   -0.3545342921352558 + 0.4753370537377221im
    0.8077682531741682 + 0.6000947508380117im
    1.0038939752615772 + 0.0im
   -0.8042569548487751 - 0.6405963628374013im
   -0.2147129522280318 + 0.9846008509750709im
   -0.2940109125982499 + 0.4580549128494589im
    0.9557804973171089 + 0.3973109143201818im
    0.9876454255450544 + 0.0im
     -0.79934335726168 - 0.6780315292599258im
   -0.0550332388392736 + 0.9592651882828983im
   -0.2760832640875706 + 0.5427315200921342im
    0.9861413269909322 + 0.14667767871633045im
    1.0506539202172076 + 0.0im
   -0.7712204915010447 - 0.6301058445223123im
   0.13552376677361944 + 0.9710647605399997im
   -0.2346148303168553 + 0.7102345469962625im
    1.0521446651366406 - 0.10416911883648566im
    0.9854255821459234 + 0.0im
   -0.7547218928194367 - 0.7307200527753428im
    0.3337597718176695 + 0.9245313060568692im
    0.7555056678669926 - 0.6921028442506652im
  -0.13962605235076483 + 0.6674441693887736im
    1.0196756885639398 - 0.403567986810368im
    1.0370385749040743 + 0.0im
   -0.7524023016981537 - 0.6554919590255408im
    0.5470631884741829 + 0.8258315816878862im
     0.856850129119688 - 0.5521031441286296im
 -0.016727201088658182 + 0.6774738851408505im
    0.8757361190397609 - 0.6977083004004283im
  -0.24985748087550927 - 0.9803418724011389im
    0.9671238563156764 + 0.0im
   -0.7711802201975106 - 0.72796856762576im
    0.7013734875241311 + 0.7038416151463572im
    0.8966012547260592 - 0.45467357068818265im
   0.11448281305898278 + 0.630316867271243im
     0.815428001504237 - 0.9930797306461538im
   -0.4735911199676696 - 0.8974524870339821im
    0.8880248242265304 + 0.0im
   -0.7775519475050343 - 0.8639629934977953im
    0.8036092545040271 + 0.5677405041111068im
    0.9516223881775341 - 0.3679142014660824im
     0.257144100072664 + 0.5984908465845797im
    0.4621345419781296 - 0.8909109795237949im
   -0.7065080130862803 - 0.7481095312530054im
    1.0307383846642757 + 0.0im
    0.9266034305156305 + 0.36637414160884907im
    0.9687174602347228 - 0.30548164125587607im
    0.4495156451294032 + 0.4715265384817793im
   -0.8712022494945793 - 0.5301357865132925im
    1.0710190218990705 + 0.0im
    -0.592733785845194 - 0.7731272168314601im
     1.009401552956874 + 0.08063846994788332im
    1.0058412648525514 - 0.25653851916514914im
    0.6878954219518458 + 0.03944322545795964im
     1.020437821976717 + 0.0im
   -0.5671526954304742 - 0.8431784218893237im
    1.0001523733486017 - 0.11991683570968216im
    0.9858810335394903 - 0.24274426776976107im
   0.47986212231230746 - 0.33375405271601283im
   -0.9933990528906721 - 0.18512496187556407im
    1.0853868963083917 + 0.0im
    0.9608752784700207 - 0.38381622577711255im
    0.9905045183602882 - 0.23075409050100096im
     0.323220770982324 - 0.5329633999824556im
   -1.0153989343358132 - 0.018062864727647357im
     1.107401832734064 + 0.0im
  -0.47976883679239213 - 0.8043568311523663im
    0.8055310647818783 - 0.6219813186081156im
    0.9870590382630294 - 0.23754586179473813im
   0.12576315738576574 - 0.5129833574641312im
   -1.0010285207387404 + 0.15063022198116155im
    1.0278846344836947 + 0.0im
   -0.4281134350046865 - 0.9132384796677624im
    0.6536686068201106 - 0.7923185512812241im
    0.9817996397568056 - 0.2622684685432363im
  -0.02726132780013805 - 0.5722711431228588im
   -0.9921227404778934 + 0.23808103174731474im
    1.0038800280720914 + 0.0im
   -0.4015537850041986 - 0.9549905976788714im
    0.4564123623164807 - 0.9057836191964054im
    0.9735371956736223 - 0.2885449655898347im
   -0.1630426221070756 - 0.7060678350480825im
   -0.9730425851918182 + 0.3078761584815582im
     0.981067803006977 + 0.0im
   -0.3290561094398114 - 0.992152384828061im
   0.25906084214213465 - 0.9578048606382372im
    0.9597267167485424 - 0.3288478275231892im
   -0.2723994522881199 - 0.6547044705098195im
    -0.963935464085029 + 0.3089163346649659im
    1.0105275761858388 + 0.0im
   -0.2638508221296871 - 0.994071343996589im
   0.15875871468724392 - 0.9054345940409332im
     0.965094992466445 - 0.33732917945443885im
  -0.34181634072887734 - 0.5865842429772671im
   -0.9604696906949206 + 0.31766529337160815im

This returns a SiteArray object which contains the gains as a flat vector with metadata about the sites, time, and frequency. To visualize the gains we can use the plotcaltable function which automatically plots the gains. Since the gains are complex we will first plot the phases.

julia
plotcaltable(angle.(intopt)) |> DisplayAs.PNG |> DisplayAs.Text

Due to the a priori calibration of the data, the gain phases are quite stable and just drift over time. Note that the one outlier around 2.5 UT is when ALMA leaves the array. As a result we shift the reference antenna to APEX on that scan causing the gain phases to jump.

We can also plot the gain amplitudes

julia
plotcaltable(abs.(intopt)) |> DisplayAs.PNG |> DisplayAs.Text

Here we find relatively stable gains for most stations. The exception is LMT which has a large offset after 2.5 UT. This is a known issue with the LMT in 2017 and is due to pointing issues. However, we can see the power of simultaneous imaging and instrument modeling as we are able to solve for the gain amplitudes and get a reasonable image.

One problem with the MAP estimate is that it does not provide uncertainties on the image. That is we are unable to statistically assess which components of the image are certain. Comrade is really a Bayesian imaging and calibration package for VLBI. Therefore, our goal is to sample from the posterior distribution of the image and instrument model. This is a very high-dimensional distribution with typically 1,000 - 100,000 parameters. To sample from this very high dimensional distribution, Comrade has an array of samplers that can be used. However, note that Comrade also satisfies the LogDensityProblems.jl interface. Therefore, you can use any package that supports LogDensityProblems.jl if you have your own fancy sampler.

For this example, we will use HMC, specifically the NUTS algorithm. For However, due to the need to sample a large number of gain parameters, constructing the posterior can take a few minutes. Therefore, for this tutorial, we will only do a quick preliminary run.

julia
using AdvancedHMC
chain = sample(rng, post, NUTS(0.8), 700; n_adapts = 500, initial_params = xopt)
PosteriorSamples
  Samples size: (700,)
  sampler used: AHMC
Mean
┌─────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────┬────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────┐
│ sky                                                                                                                                                                                                                                                                     │ instrument                                                                                                                                                                                                                                                                                                                                                                                                                         │
│ @NamedTuple{c::@NamedTuple{params::Matrix{Float64}, hyperparams::Float64}, σimg::Float64, fg::Float64}                                                                                                                                                                  │ @NamedTuple{lg::Comrade.SiteArray{Float64, 1, Vector{Float64}, Vector{Comrade.IntegrationTime{Int64, Float64}}, Vector{Comrade.FrequencyChannel{Float64, Int64}}, Vector{Symbol}}, gp::Comrade.SiteArray{Float64, 1, Vector{Float64}, Vector{Comrade.IntegrationTime{Int64, Float64}}, Vector{Comrade.FrequencyChannel{Float64, Int64}}, Vector{Symbol}}}                                                                          │
├─────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────┼────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────┤
│ (c = (params = [-0.00586159 0.00829587 … 0.00271671 -0.0253155; -0.00183311 0.000540237 … -0.0243675 -0.0213794; … ; -0.0214696 -0.0622599 … 0.0464514 0.0141755; -0.00400068 -0.0162918 … 0.0184312 0.0140476], hyperparams = 37.4628), σimg = 1.73465, fg = 0.287311) │ (lg = [0.030142, 0.0456492, 0.00531718, 0.0264892, -0.30904, 0.150679, 0.00165428, 0.0308647, -0.120561, 0.0973573  …  -0.0299363, 0.0139935, -0.729855, 0.0124434, 0.00925913, 0.0272503, -0.0974165, 0.0208916, -0.779865, 0.0110004], gp = [0.0, 0.882376, 0.0, -2.19152, 1.96554, 0.744279, 0.0, -2.24179, 2.03825, 0.556626  …  -1.33363, 0.24058, -2.18927, -1.1116, 0.0, -1.82935, -1.45448, 0.262377, -2.32636, -1.09422]) │
└─────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────┴────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────┘
Std. Dev.
┌─────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────┬──────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────┐
│ sky                                                                                                                                                                                                                                     │ instrument                                                                                                                                                                                                                                                                                                                                                                                                                   │
│ @NamedTuple{c::@NamedTuple{params::Matrix{Float64}, hyperparams::Float64}, σimg::Float64, fg::Float64}                                                                                                                                  │ @NamedTuple{lg::Comrade.SiteArray{Float64, 1, Vector{Float64}, Vector{Comrade.IntegrationTime{Int64, Float64}}, Vector{Comrade.FrequencyChannel{Float64, Int64}}, Vector{Symbol}}, gp::Comrade.SiteArray{Float64, 1, Vector{Float64}, Vector{Comrade.IntegrationTime{Int64, Float64}}, Vector{Comrade.FrequencyChannel{Float64, Int64}}, Vector{Symbol}}}                                                                    │
├─────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────┼──────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────┤
│ (c = (params = [0.54506 0.585774 … 0.547316 0.504489; 0.606577 0.656837 … 0.675332 0.563141; … ; 0.609148 0.653243 … 0.703861 0.567869; 0.584697 0.569248 … 0.611375 0.567046], hyperparams = 23.066), σimg = 0.210736, fg = 0.0862655) │ (lg = [0.145805, 0.148385, 0.0206177, 0.0222561, 0.0933108, 0.0673159, 0.0171028, 0.0188121, 0.0859491, 0.0638566  …  0.0519104, 0.0182916, 0.0969347, 0.01735, 0.0185392, 0.0208694, 0.057043, 0.0188571, 0.0989054, 0.0190878], gp = [0.0, 0.921882, 0.0, 0.017813, 0.363012, 0.934153, 0.0, 0.0172767, 0.359518, 0.94359  …  0.329351, 0.543108, 0.390215, 2.40275, 0.0, 0.0191476, 0.335076, 0.53351, 0.48877, 2.40869]) │
└─────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────┴──────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────┘

Note

The above sampler will store the samples in memory, i.e. RAM. For large models this can lead to out-of-memory issues. To avoid this we recommend using the saveto = DiskStore() kwargs which periodically saves the samples to disk limiting memory useage. You can load the chain using load_samples(diskout) where diskout is the object returned from sample.

Now we prune the adaptation phase

julia
chain = chain[(begin + 500):end]
PosteriorSamples
  Samples size: (200,)
  sampler used: AHMC
Mean
┌────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────┬────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────┐
│ sky                                                                                                                                                                                                                                                                │ instrument                                                                                                                                                                                                                                                                                                                                                                                                                             │
│ @NamedTuple{c::@NamedTuple{params::Matrix{Float64}, hyperparams::Float64}, σimg::Float64, fg::Float64}                                                                                                                                                             │ @NamedTuple{lg::Comrade.SiteArray{Float64, 1, Vector{Float64}, Vector{Comrade.IntegrationTime{Int64, Float64}}, Vector{Comrade.FrequencyChannel{Float64, Int64}}, Vector{Symbol}}, gp::Comrade.SiteArray{Float64, 1, Vector{Float64}, Vector{Comrade.IntegrationTime{Int64, Float64}}, Vector{Comrade.FrequencyChannel{Float64, Int64}}, Vector{Symbol}}}                                                                              │
├────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────┼────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────┤
│ (c = (params = [0.0516487 -0.0121713 … 0.0392888 0.0149502; -0.00316512 0.00517609 … -0.00297279 0.0287855; … ; -0.0231744 -0.0677782 … 0.0862381 0.0316802; -0.0121279 0.000865768 … 0.0474372 0.0216597], hyperparams = 38.2836), σimg = 1.78419, fg = 0.338036) │ (lg = [0.0215964, 0.0315509, 0.00658766, 0.0261082, -0.324382, 0.139737, 0.00407455, 0.0288464, -0.143645, 0.0937104  …  -0.0439703, 0.0144948, -0.783564, 0.012842, 0.00897301, 0.0280315, -0.114201, 0.0215406, -0.83272, 0.00965803], gp = [0.0, -0.336231, 0.0, -2.1916, 1.48729, -0.492485, 0.0, -2.24208, 1.56196, -0.693774  …  -1.75381, 0.470096, -2.58335, -2.55637, 0.0, -1.82805, -1.88861, 0.483017, -2.68462, -2.54654]) │
└────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────┴────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────┘
Std. Dev.
┌────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────┬───────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────┐
│ sky                                                                                                                                                                                                                                    │ instrument                                                                                                                                                                                                                                                                                                                                                                                                                        │
│ @NamedTuple{c::@NamedTuple{params::Matrix{Float64}, hyperparams::Float64}, σimg::Float64, fg::Float64}                                                                                                                                 │ @NamedTuple{lg::Comrade.SiteArray{Float64, 1, Vector{Float64}, Vector{Comrade.IntegrationTime{Int64, Float64}}, Vector{Comrade.FrequencyChannel{Float64, Int64}}, Vector{Symbol}}, gp::Comrade.SiteArray{Float64, 1, Vector{Float64}, Vector{Comrade.IntegrationTime{Int64, Float64}}, Vector{Comrade.FrequencyChannel{Float64, Int64}}, Vector{Symbol}}}                                                                         │
├────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────┼───────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────┤
│ (c = (params = [0.485001 0.557029 … 0.537258 0.508283; 0.571396 0.623497 … 0.68656 0.611878; … ; 0.635687 0.705682 … 0.764071 0.578537; 0.653908 0.629559 … 0.62431 0.56694], hyperparams = 23.9568), σimg = 0.210413, fg = 0.0784173) │ (lg = [0.152024, 0.155621, 0.0213627, 0.0220274, 0.0917437, 0.0610182, 0.0169528, 0.0190758, 0.0813932, 0.0617122  …  0.0492454, 0.020132, 0.0722768, 0.0165435, 0.0195549, 0.0223845, 0.0537558, 0.0188109, 0.0707142, 0.020175], gp = [0.0, 0.368875, 0.0, 0.0172663, 0.205655, 0.369665, 0.0, 0.0168742, 0.205041, 0.373357  …  0.191118, 0.24324, 0.44143, 0.735696, 0.0, 0.0180508, 0.195955, 0.241611, 0.735404, 0.734232]) │
└────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────┴───────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────┘

Warning

This should be run for likely an order of magnitude more steps to properly estimate expectations of the posterior

Now that we have our posterior, we can put error bars on all of our plots above. Let's start by finding the mean and standard deviation of the gain phases

julia
mchain = Comrade.rmap(mean, chain);
schain = Comrade.rmap(std, chain);

Now we can use the measurements package to automatically plot everything with error bars. First we create a caltable the same way but making sure all of our variables have errors attached to them.

julia
using Measurements
gmeas = instrumentmodel(post, (; instrument = map((x, y) -> Measurements.measurement.(x, y), mchain.instrument, schain.instrument)))
ctable_am = caltable(abs.(gmeas))
ctable_ph = caltable(angle.(gmeas))
─────────┬───────────────────┬───────────────────────────────────────────────────────────────────────────────────
      Ti │                Fr │      AA            AP          AZ          JC          LM          PV          SM
─────────┼───────────────────┼───────────────────────────────────────────────────────────────────────────────────
 0.92 hr │ 227.070703125 GHz │ 0.0±0.0       missing     missing     missing     missing  -0.34±0.37     missing
 1.22 hr │ 227.070703125 GHz │ 0.0±0.0  -2.192±0.017     missing     missing   1.49±0.21  -0.49±0.37     missing
 1.52 hr │ 227.070703125 GHz │ 0.0±0.0  -2.242±0.017     missing     missing    1.56±0.2  -0.69±0.37     missing
 1.82 hr │ 227.070703125 GHz │ 0.0±0.0  -2.276±0.018     missing     missing     1.6±0.2  -0.87±0.37     missing
 2.12 hr │ 227.070703125 GHz │ 0.0±0.0   -2.33±0.018     missing     missing    1.66±0.2  -1.05±0.37     missing
 2.45 hr │ 227.070703125 GHz │ missing       0.0±0.0     missing     missing   2.17±0.21   1.15±0.38     missing
 2.75 hr │ 227.070703125 GHz │ 0.0±0.0  -2.402±0.017     missing     missing    1.67±0.2  -1.41±0.38     missing
 3.05 hr │ 227.070703125 GHz │ 0.0±0.0  -2.432±0.019     missing     missing    1.68±0.2  -1.59±0.38     missing
 3.35 hr │ 227.070703125 GHz │ 0.0±0.0  -2.468±0.017     missing     missing   1.66±0.19  -1.78±0.38     missing
 3.68 hr │ 227.070703125 GHz │ 0.0±0.0  -2.468±0.016   1.36±0.21     missing   1.61±0.19  -1.99±0.37     missing
 3.98 hr │ 227.070703125 GHz │ 0.0±0.0  -2.439±0.019   1.22±0.21     missing   1.53±0.19  -2.16±0.51     missing
 4.28 hr │ 227.070703125 GHz │ 0.0±0.0   -2.457±0.02    1.06±0.2     missing   1.41±0.19  -2.33±0.63     missing
 4.58 hr │ 227.070703125 GHz │ 0.0±0.0  -2.373±0.019    0.88±0.2  -0.41±0.21   1.32±0.18    -2.3±1.3     missing
 4.92 hr │ 227.070703125 GHz │ 0.0±0.0  -2.423±0.021    0.68±0.2  -0.24±0.21   1.16±0.18    -1.4±2.5  -1.48±0.21
 5.18 hr │ 227.070703125 GHz │ 0.0±0.0  -2.387±0.019     0.5±0.2  -0.12±0.22   0.97±0.18   -0.41±2.8   -1.7±0.22
 5.45 hr │ 227.070703125 GHz │ 0.0±0.0  -2.305±0.022   0.34±0.19  0.023±0.22   0.75±0.18    0.22±2.9  -1.94±0.22
 5.72 hr │ 227.070703125 GHz │ 0.0±0.0       missing   0.11±0.19   0.14±0.22   0.39±0.18     missing  -2.15±0.22
 6.05 hr │ 227.070703125 GHz │ 0.0±0.0  -2.226±0.019  -0.19±0.19   0.26±0.23  -0.37±0.18     missing     missing
 6.32 hr │ 227.070703125 GHz │ 0.0±0.0  -2.164±0.019   -0.4±0.19   0.33±0.24  -1.06±0.18     missing  -2.39±0.24
 6.58 hr │ 227.070703125 GHz │ 0.0±0.0       missing  -0.68±0.19    0.4±0.24  -1.51±0.18     missing   -2.5±0.24
 6.85 hr │ 227.070703125 GHz │ 0.0±0.0  -2.109±0.021  -0.98±0.19   0.44±0.24  -1.86±0.18     missing  -2.61±0.24
 7.18 hr │ 227.070703125 GHz │ 0.0±0.0  -2.009±0.019  -1.25±0.19   0.47±0.24  -2.19±0.18     missing  -2.58±0.62
 7.45 hr │ 227.070703125 GHz │ 0.0±0.0  -1.969±0.018  -1.51±0.19   0.48±0.24  -2.41±0.18     missing  -2.55±0.84
 7.72 hr │ 227.070703125 GHz │ 0.0±0.0  -1.891±0.016  -1.75±0.19   0.47±0.24  -2.58±0.44     missing  -2.56±0.74
 7.98 hr │ 227.070703125 GHz │ 0.0±0.0  -1.828±0.018   -1.89±0.2   0.48±0.24  -2.68±0.74     missing  -2.55±0.73
─────────┴───────────────────┴───────────────────────────────────────────────────────────────────────────────────

Now let's plot the phase curves

julia
plotcaltable(ctable_ph) |> DisplayAs.PNG |> DisplayAs.Text

and now the amplitude curves

julia
plotcaltable(ctable_am) |> DisplayAs.PNG |> DisplayAs.Text

Finally let's construct some representative image reconstructions.

julia
samples = skymodel.(Ref(post), chain[begin:5:end])
imgs = intensitymap.(samples, Ref(g))

mimg = mean(imgs)
simg = std(imgs)
fig = Figure(; resolution = (700, 700));
axs = [Axis(fig[i, j], xreversed = true, aspect = 1) for i in 1:2, j in 1:2]
image!(axs[1, 1], mimg, colormap = :afmhot); axs[1, 1].title = "Mean"
image!(axs[1, 2], simg ./ (max.(mimg, 1.0e-8)), colorrange = (0.0, 2.0), colormap = :afmhot);axs[1, 2].title = "Frac. Uncer."
image!(axs[2, 1], imgs[1], colormap = :afmhot);
image!(axs[2, 2], imgs[end], colormap = :afmhot);
hidedecorations!.(axs)
fig |> DisplayAs.PNG |> DisplayAs.Text

And viola, you have just finished making a preliminary image and instrument model reconstruction. In reality, you should run the sample step for many more MCMC steps to get a reliable estimate for the reconstructed image and instrument model parameters.


This page was generated using Literate.jl.