API

Comrade

Radio Observation Sampling Exploration

abstract type AbstractCache

This defines an abstract cache that can be used to hold or precompute some computations.

abstract type AbstractModifier{M<:ComradeBase.AbstractModel} <: ComradeBase.AbstractModel

Abstract type for image modifiers. These are some model wrappers that can transform any model using simple Fourier transform properties. To see the implemented modifier

struct AddModel{T1, T2} <: Comrade.CompositeModel{T1, T2}

Adds two models together to create composite models. Note that I may change this in the future so make it easier on the compiler, i.e. make the composite model a fancy model vector and heap allocate stuff. This should help when combining multiple models together.

abstract type ArrayConfiguration

This defined the abstract type for an array configuration. Namely, baseline times, SEFD's, bandwidth, observation frequencies, etc.

abstract type CompositeModel{M1, M2} <: ComradeBase.AbstractModel

Abstract type that denotes a composite model. Where we have combined two models together.

Implementation

Any implementation of a composite type must define the following methods:

• visibility_point
• uv_combinator
• imanalytic
• visanalytic
• ComradeBase.intensity_point if model intensity is IsAnalytic
• intensitymap! if model intensity is NotAnalytic
• intensitymap if model intensity is NotAnalytic
• flux

struct ConcordanceCrescent{T} <: Comrade.GeometricModel

Creates the ConcordanceCrescent model, i.e. a flat-top crescent with a displacment and a slash. Note this creates a crescent with unit flux. If you want a different flux please use the renomed modifier.

Fields

Notes

Unlike the Gaussian and Disk models this does not create the unit version. In fact, this model could have been created using the Disk and primitives by using Comrade.jl's model composition functionality.

struct ConvolvedModel{M1, M2} <: Comrade.CompositeModel{M1, M2}

Convolves two models m1 and m2.

Notes

This is the non-exported constructor. The user should call the convolved function during use

$(TYPEDEF)

An image model given by a set of coefficients and a kernel response or basis function. This corresponds to a continous image defined by a finite set of points. The defined intensity is given by

\[ I(x,y) = \sum_{ij} c_{ij}κ(x-x_i)κ(y-y_i).\]

An important thing to note is that the $c_{ij}$ do not represent pixel intensities, i.e. the κ doesn't have to be an interpolating kernel.

Example

samples = rand(10,10)
model = DImage(samples, BSplineKernel{3})

Notes

This is defined in terms of pixel response, so the image size is 1μas. To resize the image use the scale function like with other models.

Fields

struct EHTClosurePhaseDatum{T<:Number} <: Comrade.ClosureProducts{T<:Number}

A Datum for a single closure phase. Note in the future this may get replaced with the fully covariant formalism from Blackburn et al. (2020).

struct ExtendedRing{F} <: Comrade.GeometricModel

A symmetric extended ring whose radial profile follows an inverse gamma distributions.

Note

@e mainly use this as an example of a non-analytic Fourier transform (although it has a complicated expression)

Fields

• radius

  radius of peak emission

• shape

  shape of the radial distribution

struct FFT <: Comrade.FourierTransform

Fourier transform type that specifies we will use the FFTW package to compute the Fourier transform.

Fields

• padfac

  The amount to pad the image by. Note we actually round up to the nearest factor of 2, but this will be improved in the future to use small primes

struct FFTCache{P, I} <: Comrade.AbstractCache

The cache used when the FFT algorithm is used to compute visibilties. This is an internal type and is not part of the public API

Fields

• plan

  FFTW Plan

• sitp

  FFT interpolator function

abstract type FourierTransform

Abstract type that specified which fourier transform to use

struct Gaussian{T} <: Comrade.GeometricModel

Gaussian geometrical model. This is a Gaussian with unit flux and standard deviation.

Notes

To change the Gaussian flux, and shape please use the modifier functions

abstract type GeometricModel <: ComradeBase.AbstractModel

A type that defines it is a geometric model. These are usually primitive models, and are usually analytic in Fourier and the image domain.

struct MRing{T, N} <: Comrade.GeometricModel

m-ring geometric model. This corresponds to a delta ring with a fourier expansion in θ. The m in m-ring refers to the order of the Fourier expansion. The radius is unity.

struct ModelImage{M, I, C} <: Comrade.AbstractModelImage{M}

Container for non-analytic model that contains a image cache which will hold the image, a Fourier transform cache C, which usually an instance of a <: FourierCache which usually holds a interpolator that allows you to compute visibilities.

Notes

This is an internal implementation detail that shouldn't usually be called directly. Instead the user should use the exported function modelimage, for example

using Comrade
m = ExtendedRing(20.0, 5.0)

# This creates an version where the image is dynamically specified according to the
# radial extent of the image
mimg = modelimage(m) # you can also optionally pass the number of pixels nx and ny

# Or you can create an IntensityMap
img = intensitymap(m, 100.0, 100.0, 512, 512)
mimg = modelimage(m, img)

struct NoCache <: Comrade.AbstractCache

No cache is used. This is typically used when the model is analytic in the Fourier domain.

struct ParabolicSegment{F} <: Comrade.GeometricModel

A delta parabolic segment in the image domain. The segment is centered at zero, with roots ±1 and a yintercept of 1.

ParabolicSegment(a, h)

A parabolic segment with x-intercepts ±aand a yintercept ofh``.

Note

This is just a convenience function for stretched(ParabolicSegment(), a, h)

struct PolarizedModel{TI, TQ, TU, TV} <: ComradeBase.AbstractPolarizedModel

Wrapped model for a polarized model. This uses the stokes representation of the image.

Posterior(lklhd, prior, model)

Posterior density that follows obeys the DensityInferface

The expected arguments are:

• lklhd: Which should be an intance of RadioLikelihood with whatever data products you want to fit
• prior: This should be a NamedTuple with the priors for each model ParameterHandling
• model: Function that takes a NamedTuple of parameters and constructs the Comrade model.

To evaluate the logdensity of the posterior use can either use DensityInterface.logdensityof or logdensity.

struct RenormalizedModel{M<:ComradeBase.AbstractModel, T} <: Comrade.AbstractModifier{M<:ComradeBase.AbstractModel}

Renormalizes the flux of the model to the new value flux. We have also overloaded the Base.:* operator as syntactic sugar although I may get rid of this.

struct Ring{T} <: Comrade.GeometricModel

m-ring geometric model. This corresponds to a delta ring with a fourier expansion in θ. The m in m-ring refers to the order of the Fourier expansion. The radius is unity.

struct RotatedModel{M<:ComradeBase.AbstractModel, T} <: Comrade.AbstractModifier{M<:ComradeBase.AbstractModel}

Type for the rotated model. This is more fine grained constrol of rotated model. In most use cases the end-user should be using the rotate method e.g.

rotate(model, ξ)

struct ShiftedModel{T, M<:ComradeBase.AbstractModel} <: Comrade.AbstractModifier{M<:ComradeBase.AbstractModel}

Shifts the model by Δx units in the x-direction and Δy units in the y-direction.

struct StretchedModel{M<:ComradeBase.AbstractModel, T} <: Comrade.AbstractModifier{M<:ComradeBase.AbstractModel}

Stretched the model in the x and y directions, i.e. the new intensity is

\[ I_s(x,y) = 1/(αβ) I(x/α, y/β),\]

where were renormalize the intensity to preserve the models flux.

`

struct TransformedPosterior{P<:Posterior, T}

A transformed version of the posteriorlpost`. This is an internal type that an end user shouldn't have to directly construct.

To construct a transformed posterior see the asflat, ascube, and flatten functions.

When evaluating the logdensity of the TransformedPosterior, logdenisty will include any jacobian terms automatically, and expects the argument to be the vector of transformed quantities. This is usually easier to incorporate with samplers, and optimizers, which expect homogenous vectors as arguments.

size(model)

return the size of the coefficient matrix for model.

Creates a Kamruddin and Dexter crescent model. This works by composing two disk models together.

Arguments

• router: The radius of the outer disk
• rinner: The radius of the inner disk
• shift: How much the inner disk radius is shifted (positive is to the right)
• floor: The floor of the inner disk 0 means the inner intensity is zero and 1 means it is a large disk.

_visibilities(m, u, v)

Computes the visibilities of the model m at the u,v positions u, v.

Note this is done lazily so the visibility is only computed when accessed.

amplitude(model, u, v)

Computes the visibility amplitude of model m at u,v positions u,v

If you want to compute the amplitudes at a large number of positions consider using the amplitudes function which uses MappedArrays to compute a lazy, no allocation vector.

amplitudes(m, u, v)

Computes the amplitudes of the model m at the u,v positions u, v.

Note this is done lazily so the visibility is only computed when accessed.

basemodel(model)

Returns the base model from a modified model.

bispectra(m, u1, v1, u2, v2, u3, v3)

Computes the bispectra of the model m at the triangles (u1,v1), (u2,v2), (u3,v3).

Note this is done lazily so the bispectra is only computed when accessed.

bispectrum(model, u1, v1, u2, v2, u3, v3)

Computes the complex bispectrum of model m at the uv-triangle u1,v1 -> u2,v2 -> u3,v3

If you want to compute the bispectrum over a number of triangles consider using the bispectra function which uses MappedArrays to compute a lazy, no allocation vector.

closure_phase(model, u1, v1, u2, v2, u3, v3)

Computes the closure phase of model m at the uv-triangle u1,v1 -> u2,v2 -> u3,v3

If you want to compute closure phases over a number of triangles consider using the closure_phases function which uses MappedArrays to compute a lazy, no allocation vector.

closure_phase(D1, D2, D3)

Computes the closure phase of the three visibility datums.

Notes

We currently use the high SNR Gaussian error approximation for the closure phase. In the future we may use the moment matching from Monte Carlo sampling.

closure_phases(m, u1, v1, u2, v2, u3, v3)

Computes the closure phases of the model m at the triangles (u1,v1), (u2,v2), (u3,v3).

Note this is done lazily so the closure_phases is only computed when accessed.

coherencymatrix(pimg, u, v)

Computes the coherency matrix of the polarized model.

components(m)

Returns the components for a composite model. This will return a Tuple with all the models you have constructed.

convolved(m1, m2)

Convolves two models m1 and m2. This is done lazily.

create_cache(alg, img)

Creates the model cache given for the algorithm alg using the model and a image cache image

extract_cphase(obs)

Extracts the closure phases from an ehtim observation object

Returns an EHTObservation with closure phases datums

extract_lcamp(obs)

Extracts the log-closure amp. from an ehtim observation object

Returns an EHTObservation with closure amp. datums

extract_vis(obs)

Extracts the complex visibilities from an ehtim observation object

Returns an EHTObservation with complex visibility data

load_ehtim()

Loads the eht-imaging library and stores it in the ehtim variable.

Notes

This will fail if ehtim isn't installed in the python installation that PyCall references.

load_tpy(file)

Load a ThemisPy style ascii EHT observation file.

logclosure_amplitude(model, u1, v1, u2, v2, u3, v3, u4, v4)

Computes the log-closure amplitude of model m at the uv-quadrangle u1,v1 -> u2,v2 -> u3,v3 -> u4,v3 using the formula

\[C = \log\left|\frac{V(u1,v1)V(u2,v2)}{V(u3,v3)V(u4,v4)}\right|\]

If you want to compute log closure amplitudes over a number of triangles consider using the logclosure_amplitudes function which uses MappedArrays to compute a lazy, no allocation vector.

logclosure_amplitudes(m, u1, v1, u2, v2, u3, v3, u4, v4)

Computes the log closure amplitudes of the model m at the quadrangles (u1,v1), (u2,v2), (u3,v3), (u4, v4).

Note this is done lazily so the log closure amplitude is only computed when accessed.

Construct a ModelImage where just the model m is specified

Notes

If m IsAnalytic() is the visibility domain this is a no-op and just returns the model itself. Otherwise modelimage will guess a reasonable field of view based on the radialextent function. One can optionally pass the number of pixels nx and ny in each direction.

modelimage(model, image; alg)
modelimage(m; fovx, fovy, nx, ny, pulse, alg)

Construct a ModelImage from a model, image and the optionally specified visibility algorithm alg

Notes

For analytic models this is a no-op and just return the model. For non-analytic models this wraps the model in a object with an image and precomputes the fourier transform using alg.

rad2μas(x)

Converts a number from radians to μas

renormed(model, f)

Renormalizes the model m to have total flux flux.

rotated(model, ξ)

Rotates the model by an amount ξ in radians.

shifted(model, Δx, Δy)

Shifts the model m in the image domain by an amount Δx,Δy.

smoothed(m, σ)

Smooths a model m with a Gaussian kernel with standard deviation σ.

Notes

This will created a convolved model under the hood.

stretched(model, α, β)

Stretches the model m according to the formula

\[ I_s(x,y) = 1/(αβ) I(x/α, y/β),\]

where were renormalize the intensity to preserve the models flux.

uviterator(dx, dy, nnx, nny)

Construct the u,v iterators for the Fourier transform of the image with pixel sizes dx, dy and number of pixels nx, ny

If you are extending Fourier transform stuff please use these functions to ensure that the centroid is being properly computed.

visibility(pimg, u, v)

Computes the visibility in the stokes basis of the polarized model

visibility(mimg, u, v)

Computes the complex visibility of model m at u,v positions u,v

If you want to compute the visibilities at a large number of positions consider using the visibilities function which uses MappedArrays to compute a lazy, no allocation vector.

μas2rad(x)

Converts a number from μas to rad

evpa(pimg, u, v)

electric vector position angle or EVPA of the polarized model

intensitymap!(pimg, pmodel)

Finds the polarized intensity map of the polarized model pmodel.

m̆(pimg, u, v)

Computes the fractional linear polarization in the visibility domain

\[ \breve{m} = \frac{Q + iU}{I}\]

ascube(post::Posterior)

Construct a flattened version of the posterior, where the parameters are transformed to live in the unit hypercube. In astronomy parlance, we are transforming the variables to the unit hypercube. This is done using the HypercubeTransform package.

The transformed posterior can then be evaluated by the logdensityof(transformed_posterior,x) method following the DensityInterface, where x vector that lives in the unit hypercube. Note this already includes the jacobian of the transformation so this does not need to be added.

This transform is useful for NestedSampling methods that often assume that the model is written to live in the unit hypercube.

asflat(post::Posterior)

Construct a flattened version of the posterior, where the parameters are transformed so that their support is from (-∞, ∞). This uses TransformVariables

The transformed posterior can then be evaluated by the logdensityof(transformed_posterior,x) method following the DensityInterface, where x is a flattened vector of the infinite support variables. Note this already includes the jacobian of the transformation so this does not need to be added.

This is useful for optimization and sampling algorithms such as HMC that will use gradients to explore the posterior surface.

flatten(post::Posterior)

Flatten the representation of the posterior. Internally this uses ParameterHandling to construct a flattened version of the posterior.

Note this is distinct from asflat that transforms the variables to live in (-∞,∞). Instead this method just flattens the repsentation of the model from a NamedTuple to a vector. This allows the easier integration to optimization and sampling algorithms.

plot(image::IntensityMap)

where image is templated off of EHTImage struct.

Details

This was created to be close to the ehtim display object. It takes an EHTImage object and plots it according to EHT conventions.

Note that is does not save the figure.

plot(image::IntensityMap)

where image is templated off of EHTImage struct.

Details

This was created to be close to the ehtim display object. It takes an EHTImage object and plots it according to EHT conventions.

Note that is does not save the figure.

inverse(posterior::TransformedPosterior, x)

Transforms the value model parameters x into the flattened transformed space.

transform(posterior::TransformedPosterior, x)

Transforms the value x into parameter space, i.e. usually a NamedTuple.