API

DiskWindModel(r̄::Float64, rFac::Float64, α::Float64, i::Float64; rot::Float64=0.0, 
        nr::Int=128, nϕ::Int=256, scale::Symbol=:log, kwargs...)

Uses the model constructor to create a DiskWind model of the BLR as detailed in Long+2023 and Long+2025.

Parameters

• r̄::Float64: Mean radius of model (in terms of $r_s$)
• rFac::Float64: Radius factor
• α::Float64: Power-law source function scaling
• i::Float64: Inclination angle in radians
• rot::Float64=0.0: Rotation of system plane about z-axis in radians
• nr::Int=128: Number of radial bins
• nϕ::Int=256: Number of azimuthal bins
• scale::Symbol=:log: Radial binning scale (:log or :linear)
• kwargs...: Extra keyword arguments for model constructor (see examples)

Returns

• model object

Note

Similar to another DiskWind model constructor but here we pass r̄, rFac, and α.

DiskWindModel(rMin::Float64, rMax::Float64, i::Float64; nr::Int=128, nϕ::Int=256, 
        I::Function=DiskWindIntensity, v::Function=vCircularDisk, scale::Symbol=:log, kwargs...)

Uses the model constructor to create a DiskWind model of the BLR as detailed in Long+2023 and Long+2025.

Parameters

• rMin::Float64: Minimum radius of model (in terms of $r_s$)
• rMax::Float64: Maximum radius of model (in terms of $r_s$)
• i::Float64: Inclination angle in radians (all rings have the same inclination)
• nr::Int=128: Number of radial bins
• nϕ::Int=256: Number of azimuthal bins
• I::Function=DiskWindIntensity: Intensity function
• v::Function=vCircularDisk: Velocity function
• scale::Symbol=:log: Radial binning scale (:log or :linear)
• kwargs...: Extra keyword arguments for I and v functions (see examples)

Returns

• model object

Note

Similar to other DiskWind model constructor but must explicitly pass rMin and rMax.

cloudModel(nClouds::Int64; μ::Float64=500., β::Float64=1.0, F::Float64=0.5, rₛ::Float64=1.0, θₒ::Float64=π/2, γ::Float64=1.0, ξ::Float64=1.0, i::Float64=0.0, I::Union{Function,Float64}=IsotropicIntensity, v::Union{Function,Float64}=vCircularCloud, rng::Union{AbstractRNG,Symbol}=Random.GLOBAL_RNG, seed::Union{Nothing,Integer}=nothing, parallel::Bool=true, kwargs...)

Uses the model constructor to create a cloud model of the BLR similar to Pancoast+ 2011 and 2014.

Parameters

• nClouds::Int64: Number of clouds
• μ::Float64=500.: Mean radius of model (in terms of $r_s$)
• β::Float64=1.0: Shape parameter for radial distribution
• F::Float64=0.5: Beginning radius in units of μ where clouds can be placed
• rₛ::Float64=1.0: Scale radius (in terms of $r_s$)
• θₒ::Float64=π/2: Maximum opening angle of cloud distribution (rad)
• γ::Float64=1.0: Disk concentration parameter
• ξ::Float64=1.0: Fraction of clouds in back side that have not been moved to the front (when ξ = 1.0 clouds equally distributed front - back and when ξ = 0.0 all clouds are on the front side)
• i::Float64=0.0: Inclination angle of system (rad)
• I::Union{Function,Float64}=IsotropicIntensity: Intensity function
• v::Union{Function,Float64}=vCircularCloud: Velocity function
• rng::Union{AbstractRNG,Symbol}=Random.GLOBAL_RNG: Random number generator. The default AbstractRNG path preserves the legacy sequential draw order; rng=:philox uses independent counter-based per-cloud streams.
• seed::Union{Nothing,Integer}=nothing: Required when rng=:philox. For legacy seeded models, pass rng=MersenneTwister(seed).
• parallel::Bool=true: Whether to use threaded cloud generation for the rng=:philox path.
• kwargs...: Extra keyword arguments for I and v functions (see examples)

Returns

• model object

Note

Similar to the other cloudModel methods, but here the per-cloud ϕ₀, rot, and θₒ values are drawn internally (one independent stream per cloud) while i is held constant for the system.

cloudModel(ϕ₀::Vector{Float64}, i::Vector{Float64}, rot::Vector{Float64}, θₒ::Vector{Float64}, 
        θₒSystem::Float64, ξ::Float64; rₛ::Float64=1.0, μ::Float64=500., β::Float64=1.0, F::Float64=0.5, 
        I::Union{Function,Float64}=IsotropicIntensity, v::Union{Function,Float64}=vCircularCloud, kwargs...)

Uses the model constructor to create a cloud model of the BLR similar to Pancoast+ 2011 and 2014.

Parameters

• ϕ₀::Vector{Float64}: Initial azimuthal angle of cloud (rad)
• i::Vector{Float64}: Inclination angle (rad)
• rot::Vector{Float64}: Random rotation of cloud about z axis (rad)
• θₒ::Vector{Float64}: Opening angle of cloud (rad)
• θₒSystem::Float64: Maximum opening angle of the system (rad)
• ξ::Float64: Fraction of clouds in back side that have not been moved to the front (when ξ = 1.0 clouds equally distributed front - back and when ξ = 0.0 all clouds are on the front side)
• rₛ::Float64=1.0: Scale radius (in terms of $r_s$)
• μ::Float64=500.: Mean radius of model (in terms of $r_s$)
• β::Float64=1.0: Shape parameter for radial distribution
• F::Float64=0.5: Beginning radius in units of μ where clouds can be placed.
• I::Union{Function,Float64}=IsotropicIntensity: Intensity function
• v::Union{Function,Float64}=vCircularCloud: Velocity function
• kwargs...: Extra keyword arguments for I and v functions (see examples)

Returns

• model object

Note

Similar to other cloudModel method but here you must explicitly pass ϕ₀, i, rot, and θₒ.

fillDiskGrid!(rSystem, ϕSystem, ϕ₀, η, xSystem, ySystem, zSystem, α, β,
    i, rot, θₒ, M, r3D, rMin, rMax, ηₒ, η₁, αRM, rNorm)

Typed inner pixel loop of the disk-wind model constructor (function barrier).

The constructor body is type-unstable (its locals come from kwargs), which would box every per-pixel value; routing the loop through this concretely-typed helper keeps the raytracing allocation-free. Fills the supplied matrices in place; out-of-range pixels (r < rMin or r > rMax, compared after rounding to 9 significant digits) are set to NaN. M and r3D are the precomputed matrices described in raytrace; ηₒ, η₁, αRM, rNorm are the response parameters (hoisted from kwargs once instead of splatted per pixel).

Rows are processed in parallel when Julia is started with multiple threads (julia -t N). Every pixel writes only its own matrix slots and the computation is deterministic, so results are bit-identical at any thread count.

camera

Camera coordinates struct.

Fields

• α::Union{Vector{Float64}, Matrix{Float64}}: x values in the camera plane
• β::Union{Vector{Float64}, Matrix{Float64}}: y values in the camera plane
• raytraced::Bool: whether the camera has been used to raytrace the model

model

A mutable structure to hold many rings and their parameters that model the BLR.

Fields

• rings::Vector{ring}: List of ring objects, see ring struct
• profiles::Union{Nothing,Dict{Symbol,profile}}: Dictionary of profiles (see profile struct) with keys as symbols; optional, usually initialized to empty dictionary and filled in with setProfile!
• camera::Union{Nothing,camera}: Camera coordinates (α,β) corresponding to each ring used to generate images and in raytracing, see camera struct
• subModelStartInds::Vector{Int}: Indices of start of each submodel in list of rings; used to separate out submodels for raytracing or for the recovery of individual models after being combined
• cache::Union{Nothing,Dict{Any,Array}}: Memoized getVariable results so repeated profile calculations don't re-gather data from the rings. Managed automatically by the package's own mutating functions; if you mutate ring fields directly, call reset!(m) afterwards to invalidate it. Set to nothing to disable caching for a model entirely.
• params::Union{Nothing,NamedTuple}: Construction provenance recorded by the public constructors (DiskWindModel, cloudModel) and propagated recursively through + and raytrace!. Its first entry is always constructor=<Symbol> naming how the model was built; leaf records store the arguments that constructor was called with, :+ records left/right subtrees, and :raytrace! records the raytrace arguments plus a parent subtree. Used by rebuild to reconstruct (or re-parameterize) a model. nothing when the model was not built by a public constructor. This captures construction history ONLY – in-place mutators (removeNaN!, zeroDiskObscuredClouds!) are NOT recorded, so a caller who ran them must re-apply the same calls in the same order to a rebuilt model. removeNaN! is observable-neutral (every binned observable already skips non-finite points), so the only unrecorded mutator with a physical effect is zeroDiskObscuredClouds!, which is deterministic and cheap to re-run. params is metadata only and never affects numerics.

Constructors

model(rings::Vector{ring}, profiles::Union{Nothing,Dict{Symbol,profile}}, camera::Union{Nothing,camera}, subModelStartInds::Vector{Int}) 

• the most flexible constructor, takes in user supplied vector of ring objects, optional dictionary of profiles, optional camera coordinates, and vector of indices for submodels

model(rings::Vector{ring{Vector{Float64},Float64}})

• cloud model constructor, takes in a vector of ring objects that are each points in space and returns a model object with camera coordinates calculated from the physical parameters of each ring

model(rMin::Float64, rMax::Float64, i::Float64, nr::Int, nϕ::Int, I::Function, v::Function, scale::Symbol; kwargs...)

• disk-wind model constructor, takes in minimum and maximum radius, inclination angle, number of radial and azimuthal bins, intensity and velocity functions, and scale for radial binning; returns a model object with rings generated from these parameters

model(r̄::Float64, rFac::Float64, Sα::Float64, i::Float64, nr::Int, nϕ::Int, scale::Symbol; kwargs...)

• disk-wind model constructor, takes in average radius, radius scaling factor, power law for source function, inclination angle, number of radial and azimuthal bins, and scale for radial binning; returns a model object with rings generated from these parameters

profile

A struct to hold binned data, usually bound to model struct with profiles.jl#setProfile!.

Fields

• name::Symbol: Name of profile
• binCenters::Vector{Float64}: Bin centers
• binEdges::Vector{Float64}: Bin edges
• binSums::Vector{Float64}: Sum of values in each bin (crude integral over bin)

ring{V,F} <: AbstractRing{V,F}

A mutable structure to hold parameters of each model ring, where the "ring" is a circle in the camera plane observing the BLR.

Fields

• r: Distance from central mass (in terms of $r_s$)

  • Union{Vector{Float64}, Float64, Function}
  • Can be a single value for constant radius, or vector corresponding to azimuthal angles
  • Can be a function returning Vector{Float64} or Float64

• i: Inclination angle in radians

  • Union{Vector{Float64}, Float64}
  • Must be between 0 and $\pi/2$, with 0 being face-on and $\pi/2$ being edge-on

• rot: Rotation of system plane about z-axis in radians

  • Union{Vector{Float64}, Float64}

• θₒ: Opening angle of ring in radians

  • Union{Vector{Float64}, Float64}
  • Should be between 0 and $\pi/2$

• v: Line of sight velocity

  • Union{Vector{Float64}, Float64, Function}
  • Can be a function that calculates velocity from other parameters

• I: Intensity

  • Union{Vector{Float64}, Float64, Matrix{Float64}, Function}
  • Can be a function that calculates intensity from other parameters

• ϕ: Azimuthal angle in radians

  • Union{Vector{Float64}, Float64}

• ϕ₀: Initial azimuthal angle before rotation in radians

  • Union{Vector{Float64}, Float64}
  • Defaults to 0.0 if not provided or if rot is 0.0

• ΔA: Projected area of each ring element in image

  • Union{Vector{Float64}, Float64}
  • Used in calculating profiles

• reflect: Whether cloud is reflected across disk mid-plane

  • Union{Bool, Array{Bool,}}

• τ: Optical depth

  • Union{Vector{Float64}, Float64, Function}

• η: Response parameter for reverberation

  • Union{Vector{Float64}, Float64, Function}

• Δr: Distance between camera pixels in r

  • Float64

• Δϕ: Distance between camera pixels in ϕ

  • Float64

• scale: Encoding for camera ring scaling

  • Union{Nothing, Symbol}
  • :log or :linear scale

• x, y, z: Cached 3D system coordinates of each point (camera at +x), or nothing

  • Union{Vector{Float64}, Float64, Nothing}
  • Computed once at construction (or lazily by getXYZ) so hot loops never recompute rotate3D
  • Values are post-reflection: if reflect is true the stored coordinates already include it
  • Default nothing; access through getXYZ(ring) which fills the cache on first use

Constructor

ring(; r, i, v, I, ϕ, rot=0.0, θₒ=0.0, ϕ₀=0.0, ΔA=1.0, reflect=false,
    τ=0.0, η=1.0, Δr=1.0, Δϕ=1.0, scale=nothing, x=nothing, y=nothing, z=nothing, kwargs...)

Required parameters:

• r, i, v, I, ϕ

Optional parameters with defaults:

• rot=0.0, θₒ=0.0, ϕ₀=0.0, ΔA=1.0, reflect=false, τ=0.0, η=1.0, Δr=1.0, Δϕ=1.0, scale=nothing

Additional keyword arguments are passed to the r, v, I, and τ functions if they are provided as functions.

drawCloud(;μ::Float64=500.,β::Float64=1.0,F::Float64=0.5,
        ϕ₀::Float64=0.0,i::Float64=π/4,rot::Float64=0.0,rₛ::Float64=1.0,
        θₒ::Float64=0.0,θₒSystem::Float64=0.0,I::Union{Function,Float64}=IsotropicIntensity,
        v::Union{Function,Float64}=vCircularCloud,ξ::Float64=1.0,
        rng::AbstractRNG=Random.GLOBAL_RNG,kwargs...)

Generates a model ring struct for a single cloud drawn from a thick disk-like structure with parameters defined by the input arguments similar to Pancoast+ 2011 and 2014.

Arguments:

• μ::Float64=500.: Mean radius of model (in terms of $r_s$)
• β::Float64=1.0: Shape parameter for radial distribution
• F::Float64=0.5: Beginning radius in units of μ where clouds can be placed.
• ϕ₀::Float64=0.0: Initial azimuthal angle of cloud (rad)
• i::Float64=π/4: Inclination angle (rad)
• rot::Float64=0.0: Random rotation of cloud about z axis (rad)
• rₛ::Float64=1.0: Scale radius (in terms of $r_s$)
• θₒ::Float64=0.0: Opening angle of cloud (rad)
• θₒSystem::Float64=0.0: Maximum opening angle of the system (rad)
• ξ::Float64=1.0: Fraction of clouds in back side that have not been moved to the front (when ξ = 1.0 clouds equally distributed front - back and when ξ = 0.0 all clouds are on the front side)
• I::Union{Function,Float64}=IsotropicIntensity: Intensity function/value
• v::Union{Function,Float64}=vCircularCloud: Velocity function/value
• kwargs...: Extra keyword arguments for I and v functions (see examples)

Returns:

A model ring struct representing the properties of the cloud.

getG(β::Float64)

Returns a Gamma distribution with parameters derived from the input parameter β as in Pancoast+ 2014 equation 12.

getGamma(;μ::Float64,β::Float64,F::Float64)

Returns a Gamma distribution with parameters derived from the input parameters as in Pancoast+ 2014 equations 7-10.

getR(rₛ::Float64,μ::Float64,β::Float64,F::Float64,
    g::Gamma{Float64},rng::AbstractRNG=Random.GLOBAL_RNG)

Returns the radius r calculated using Pancoast+ 2014 equation 12, where g is a Gamma distribution and a random number is drawn from it.

getR(rₛ::Float64,γ::Gamma{Float64},rng::AbstractRNG=Random.GLOBAL_RNG)

Returns a random number drawn from the Gamma distribution γ and shifted by rₛ.

getR(rₛ::Float64,μ::Float64,β::Float64,F::Float64,g::Float64)

Returns the radius r calculated using Pancoast+ 2014 equation 12, where g is already set.

Base.:+(m1::model, m2::model) combines two models by concatenating their rings and camera parameters.

Create a new model with the combined rings and camera parameters, and updates the subModelStartInds accordingly.

Base.getindex(m::model, i::Int)

Retrieves the i-th submodel from the model m.

The extracted submodel carries the matching side of a :+ params provenance tree when that side can be identified unambiguously (so e.g. (a+b)[2].params == b.params and the result is rebuild-able); otherwise its params is nothing – in particular for operands built without a public constructor and for submodels sliced out of a multi-slot raytrace!d model (re-raytracing a slice alone is not equivalent to slicing the raytraced whole, so no faithful construction record exists).

Base.getindex(cm::CompositeModel, line::String) -> model

Retrieve the model registered for line. Errors (listing the known line names) if line is not present.

Base.length(cm::CompositeModel) -> Int

Number of lines registered in cm.

Base.show(io::IO, cm::CompositeModel)

Print one row per line: name, lineCenter, fluxRatio, number of rings, and whether profiles have been computed (setProfile!/getProfile) for that line's model.

_finiteVRange(m::model) -> (vmin, vmax)

Internal helper: the definition of a line's velocity range – the min/max of v over points where isfinite(v) && isfinite(I*ΔA) (the W4-T4 finiteness rule; finite-I alone would let a NaN v poison the range). Single non-allocating pass over the memoized flattened v/I/ΔA gathers. Shared by lineOverlap, getSpectrum (and Workstream 5) so the three cannot disagree.

Returns (NaN, NaN) when the model has no finite points; callers treat that as "this line contributes no range" (skip it).

_fluxWeights(cm::CompositeModel) -> Dict{String,Float64}

Internal helper: per-line multiplicative weight fluxRatios[line] / totalFlux(line) that renormalizes each line's binned profile to unit integral before scaling by its fluxRatio – the integrated-flux semantics of decision 3. totalFlux(line) = Σ I[k]*ΔA[k] over points where isfinite(v[k]) && isfinite(I[k]*ΔA[k]) – both conditions, matching exactly what binAccumulate! counts (finite-I/NaN-position points occur in practice, so finite-I alone is not enough).

Uses the memoized flattened gathers (getVariable(m, …; flatten=true) hits model.cache) and one non-allocating loop. Deliberately not memoized itself – it recomputes per call so nothing new joins the model.cache invalidation contract; the O(N) sum is noise next to the binning it accompanies.

Errors (naming the offending line) if any line's total flux is not positive: a zero-emission line cannot be normalized to unit integral, and the alternative – an Inf weight whose non-finite weighted fluxes the binning silently drops – would make the line quietly vanish from the spectrum.

addLine!(cm::CompositeModel, m::model; line::String, lineCenter::Float64, fluxRatio::Float64=1.0) -> cm

Register an independently built model (which may itself be a +-combined, raytraced multi-component model) as a new line in cm. Returns cm.

Parameters

• line::String: unique (case-sensitive) name for this line.
• lineCenter::Float64: central wavelength (same units convention as cm's other lines). Must be positive.
• fluxRatio::Float64=1.0: this line's velocity-integrated flux relative to cm.lines[1]. Must be positive.

addLine!(cm::CompositeModel; line::String, lineCenter::Float64, fluxRatio::Float64=1.0,
         from::String=cm.lines[1], overrides...) -> cm

Add a new line by parameter reuse: rebuild the recorded construction parameters (model.params, see the model docstring and Workstream 4-T1) of the from line, merged with overrides, and register the result under line. Because rebuild handles nested :+/:raytrace! records, this works for combined and raytraced from lines too – overrides broadcast to every leaf submodel that has a matching parameter (so e.g. overriding τ on a raytraced disk+cloud line reaches both submodels' intensity functions, if both accept τ).

Errors with a clear message if the from model has no recorded params (params === nothing, i.e. it was not built by a public constructor) – pass an explicit model via the other addLine! method instead.

Override semantics – constructor semantics, no geometry invariance

rebuild does exactly what calling the original public constructor with the merged arguments would do, nothing more – no geometry invariance is promised. Whether an override moves geometry depends on which parameterization was recorded:

• :DiskWindModelMean (i.e. DiskWindModel(r̄, rFac, α, i; ...)): α (and r̄, rFac) set rMin/rMax via get_rMinMaxDiskWind, so overriding any of them moves the whole grid.
• :DiskWindModel (i.e. explicit DiskWindModel(rMin, rMax, i; ...)): α reaches only the intensity function via kwargs – geometry (r, ϕ, v) is untouched, only I changes.

Cloud-model seed rule

The stored rng/seed live at the leaves of a (possibly nested) params tree, so this rule is applied per leaf, recursively:

• rng=:philox and the caller does not override seed -> the stored seed is reused automatically (unmentioned override keys are left as recorded), producing identical clouds.
• rng=:philox and the caller overrides seed -> fresh gas from the new seed. Overriding with seed=nothing is rejected with a clear error (cloudModel(...; rng=:philox) requires an explicit integer seed) – philox cannot draw without one.
• rng=:legacy (the model consumed Random.GLOBAL_RNG at construction time) -> a single @warn is emitted noting that gas geometry cannot be reproduced exactly and will be re-drawn (statistically equivalent, not identical) from GLOBAL_RNG.

Not recorded

params captures construction history only. In-place mutators (removeNaN!, zeroDiskObscuredClouds!) are not recorded and are not re-applied by rebuild/ addLine! – see the rebuild docstring's "Not recorded" section. A caller who ran such mutators on the from model must re-apply the same calls, in the same order, to the newly added line's model.

Parameters

• line::String: unique (case-sensitive) name for the new line.
• lineCenter::Float64: central wavelength. Must be positive.
• fluxRatio::Float64=1.0: this line's flux relative to cm.lines[1]. Must be positive.
• from::String=cm.lines[1]: which registered line's params to reuse.
• overrides...: forwarded to rebuild (constructor-argument overrides).

getSpectrum(cm::CompositeModel; bins=100, z=0.0, kwargs...)
    -> (edges, centers, flux::Dict{String,Vector{Float64}}, total::Vector{Float64})

Combined wavelength-space spectrum of all lines. Each line's flattened points are binned by wavelength(v, lineCenter)*(1+z) weighted by I*ΔA*_fluxWeights(cm)[line], over a shared wavelength grid spanning every line's finite range (from _finiteVRange) times (1+z). total is the elementwise sum of the per-line vectors. With the unit-integral weighting, sum(flux[line]) equals fluxRatios[line] (decision 3).

Arguments

• bins::Union{Int,Vector{Float64}}=100: bin count, or an explicit edge vector (passthrough as in binnedSum). For an integer count no edge vector is materialized – the same minX/maxX go to every line's binnedSum, whose deterministic constructBinEdges yields bit-identical edges and preserves the uniform direct-index fast path.
• z::Real=0.0: redshift; z=0 = rest frame, z>0 shifts to the observed frame.
• minX/maxX::Union{Real,Nothing}=nothing: pin one or both ends of the shared wavelength grid (observed frame, i.e. after the (1+z) shift) instead of auto-computing them from the lines' finite ranges – forwarded to every line's binnedSum like the other binning keywords the simple model methods forward. Ignored when bins is an edge vector (the edges pin the grid).
• centered::Bool=false: forwarded to binnedSum/constructBinEdges. The default false builds edges at exactly [minX, maxX]; centered=true pads them half a bin so the extremes fall in bin centers, as in the getProfile default.
• overflow::Union{Bool,Nothing}=nothing: an explicit true/false is forwarded to binnedSum untouched. The default (nothing) resolves to true for integer bins and false for a user-supplied edge vector: non-centered internal edges built at exactly [λmin, λmax] place each line's extremal points on the boundary, which binnedSum otherwise drops (deliberate package-wide convention) – overflow=true is the sanctioned mechanism to keep them, is what makes the integrated flux exact, and (for a user-narrowed minX/maxX window) keeps out-of-window flux by accumulating it into the boundary bins. Pass overflow=false explicitly for drop-outside-the-window semantics; against a user edge vector (a data grid) out-of-range flux drops by default.

getVariable(cm::CompositeModel, variable; line::String, flatten=false)

Per-line forwarding of getVariable: gather variable from the model registered under line.

gpu(cm::CompositeModel; T=Float32, kwargs...) -> ResidentCompositeModel

Move every line's model onto the GPU with the per-model gpu (each line gets the same keyword arguments, e.g. T=Float32) and wrap the per-line ResidentModel handles as a ResidentCompositeModel. Requires CUDA.jl to be loaded so the package extension can provide the CUDA-backed gpu(::model) method — without it the gpu(::Any) fallback error fires per line.

image(cm::CompositeModel, [variable]; line::String, kwargs...)

Per-line forwarding of the image recipe: renders cm[line] exactly as image(cm[line], variable; kwargs...) would. line is required (no default) – unlike plot3d, an intensity image has no meaningful "overlay all lines" mode (each line's I is in its own arbitrary units, so mixing them in one color scale would be misleading).

image(rcm::ResidentCompositeModel, [variable]; line::String, kwargs...)

Per-line forwarding of the image recipe for a ResidentCompositeModel: renders rcm[line] exactly as image(rcm[line], variable; kwargs...) would (the ResidentModel recipe path – host copies of that line's device columns). As on the host composite, line is required (no default): an intensity image has no meaningful "overlay all lines" mode. variable must be a ModelArrays column name (custom Function variables need the host model path).

lineOverlap(cm::CompositeModel) -> Vector{NamedTuple}

Report wavelength-space overlaps between lines. Each line spans [lineCenter*(1+vmin), lineCenter*(1+vmax)] from _finiteVRange; every pair whose intervals intersect yields an entry (lineA, lineB, λlo, λhi) (the overlap interval [λlo, λhi]). An empty vector means no lines overlap. Lines with no finite points (_finiteVRange = NaN) are skipped.

A ResidentCompositeModel is supported through the same shared implementation (the ::ResidentCompositeModel method lives with its resident siblings in gpu_observables.jl, W4-G3); the only difference is that each line's velocity range then comes from the device reductions of _finiteVRange(::ResidentModel) instead of the host gathers.

lineRatio(cm::CompositeModel, a::String, b::String) -> Float64

Integrated flux ratio of line a to line b. By the integrated-flux semantics of fluxRatio (each line's profile is normalized to unit integral before scaling, see _fluxWeights) this is exactly cm.fluxRatios[a] / cm.fluxRatios[b] – the constant that the velocity-resolved getProfile(cm, :ratio; lines=(a, b)) profile integrates to (and equals bin-by-bin in the degenerate identical-geometry case). Errors (listing the known lines) if either name is not registered.

A ResidentCompositeModel is supported through the same shared implementation (the ::ResidentCompositeModel method lives with its resident siblings in gpu_observables.jl, W5-G1) – fluxRatios is host metadata on both, so the two methods are identical.

plot3d(cm::CompositeModel, [variable], [annotate]; line=nothing, kwargs...)

Per-line forwarding of the plot3d recipe for a CompositeModel.

• line::String: identical to plot3d(cm[line], variable, annotate; kwargs...).
• line=nothing (default): overlay every line's geometry on one 3D plot, one solid color per line. Reuses _plot3d_submodels – the same per-submodel point gather the single-model plot3d recipe uses – for each line's (x,y,z) (and its NaN mask on variable); unlike the single-model recipe this does NOT color points by variable's value (a shared colormap across lines with independently-scaled intensities/variables would not be meaningful) – points are colored by line instead, and variable only selects which column supplies the per-point finite/NaN mask (default geometry, i.e. no masking, matching the single-model recipe's default).

plot3d(rcm::ResidentCompositeModel, [variable], [annotate]; line=nothing, kwargs...)

Per-line forwarding of the plot3d recipe for a ResidentCompositeModel, matching plot3d(::CompositeModel, …):

• line::String: identical to plot3d(rcm[line], variable, annotate; kwargs...) (the ResidentModel recipe path).
• line=nothing (default): overlay every line's geometry on one 3D plot, one solid color per line, via the same shared implementation as the host composite overlay. The _plot3d_submodels(::ResidentModel, …) restrictions apply per line: a multi-submodel line must carry raytrace metadata (build the handle with raytrace=true), and custom-Function mask variables are not supported (pass a column Symbol, or use the host model).

raytrace!(cm::CompositeModel; line=nothing, kwargs...) -> cm

Per-line forwarding of raytrace!. line=nothing (the default) raytraces every line; otherwise only the named line.

raytrace! returns a new model (it does not mutate its argument), so this reassigns cm.models[line] = raytrace!(cm.models[line]; kwargs...). Lines that raytrace! would warn about and hand back unaltered – a single-submodel model (subModelStartInds == [1]) or an already-raytraced one (camera.raytraced) – are skipped silently so a multi-line loop does not warn-spam.

Raytracing happens within each line's own submodels only. Cross-line occlusion is deliberately not modeled: different lines are at different wavelengths and each line's τ is that line's own optical depth (see the Workstream-4 plan, "Out of scope").

removeNaN!(cm::CompositeModel) -> cm

Per-line forwarding of removeNaN!: drop NaN-sentinel points from every line's model.

reset!(cm::CompositeModel; line=nothing, profiles=true, img=false) -> cm

Per-line forwarding of reset!. line=nothing (the default) resets every line; otherwise only the named line. Extra keywords (profiles, img) forward to the per-line model method.

resident(cm::CompositeModel; T=Float64, backend=KernelAbstractions.CPU(), raytrace=false)
    -> ResidentCompositeModel

Flatten every line's model with the per-model resident (each line gets the same keyword arguments) and wrap the per-line ResidentModel handles as a ResidentCompositeModel. The default CPU() backend keeps everything on the host — useful for testing the resident composite pipeline without a GPU. Use gpu(cm) (with CUDA.jl loaded) to build a device-resident composite.

setProfile!(cm::CompositeModel, p::profile; line::String, overwrite=false) -> cm

Per-line forwarding of setProfile!: store profile p on the model registered under line.

spectrum(cm::CompositeModel; z=0.0, bins=100, kwargs...)
spectrum(rcm::ResidentCompositeModel; z=0.0, bins=100, kwargs...)

Combined wavelength-space spectrum plot: one series per line plus a black "total" series (both from getSpectrum), with pairwise line overlaps (lineOverlap) shaded as gray bands behind the curves.

A ResidentCompositeModel works through the same recipe body with no extra branch (W4-G3): getSpectrum and lineOverlap resolve at call time and dispatch to the resident methods (gpu_observables.jl), which run the per-line reductions/binning on the device. The resident getSpectrum restrictions apply (a user-supplied bins edge vector must be uniform).

Keywords

• z::Real=0.0: redshift, forwarded to getSpectrum (z=0 = rest frame; matches getSpectrum's observed-frame convention).
• bins::Union{Int,Vector{Float64}}=100: bin count or explicit edge vector, forwarded to getSpectrum.
• other kwargs...: ordinary Plots attributes.

lineOverlap reports rest-frame wavelength intervals; they are scaled by (1+z) here before shading so the bands line up with the (possibly redshifted) spectrum.

spectrum(cm::CompositeModel; z=0.0, bins=100, kwargs...)
spectrum(rcm::ResidentCompositeModel; z=0.0, bins=100, kwargs...)

Combined wavelength-space spectrum plot: one series per line plus a black "total" series (both from getSpectrum), with pairwise line overlaps (lineOverlap) shaded as gray bands behind the curves.

A ResidentCompositeModel works through the same recipe body with no extra branch (W4-G3): getSpectrum and lineOverlap resolve at call time and dispatch to the resident methods (gpu_observables.jl), which run the per-line reductions/binning on the device. The resident getSpectrum restrictions apply (a user-supplied bins edge vector must be uniform).

Keywords

• z::Real=0.0: redshift, forwarded to getSpectrum (z=0 = rest frame; matches getSpectrum's observed-frame convention).
• bins::Union{Int,Vector{Float64}}=100: bin count or explicit edge vector, forwarded to getSpectrum.
• other kwargs...: ordinary Plots attributes.

lineOverlap reports rest-frame wavelength intervals; they are scaled by (1+z) here before shading so the bands line up with the (possibly redshifted) spectrum.

spectrum(cm::CompositeModel; z=0.0, bins=100, kwargs...)
spectrum(rcm::ResidentCompositeModel; z=0.0, bins=100, kwargs...)

Combined wavelength-space spectrum plot: one series per line plus a black "total" series (both from getSpectrum), with pairwise line overlaps (lineOverlap) shaded as gray bands behind the curves.

A ResidentCompositeModel works through the same recipe body with no extra branch (W4-G3): getSpectrum and lineOverlap resolve at call time and dispatch to the resident methods (gpu_observables.jl), which run the per-line reductions/binning on the device. The resident getSpectrum restrictions apply (a user-supplied bins edge vector must be uniform).

Keywords

• z::Real=0.0: redshift, forwarded to getSpectrum (z=0 = rest frame; matches getSpectrum's observed-frame convention).
• bins::Union{Int,Vector{Float64}}=100: bin count or explicit edge vector, forwarded to getSpectrum.
• other kwargs...: ordinary Plots attributes.

lineOverlap reports rest-frame wavelength intervals; they are scaled by (1+z) here before shading so the bands line up with the (possibly redshifted) spectrum.

wavelength(v, lineCenter) -> λ

First-order velocity→wavelength map λ = lineCenter * (1 + v) (decision 4), with v the stored line-of-sight velocity in units of c. Works for scalar or vector v (broadcasts). Stored v is redshift-positive (positive v = receding), so (1 + v) is the correct sign – do not re-derive it.

CompositeModel

Holds an ordered collection of line => model pairs so that several broad emission lines (e.g. Hα + Hβ + CIV) can be modeled together. The existing model/ring/ raytrace!/cache/GPU code paths are unmodified and know nothing about lines – CompositeModel is purely a wrapper around independently valid model objects.

Fields

• lines::Vector{String}: unique line names, in creation order. lines[1] is the reference line (the implicit from default for addLine! parameter reuse, and the line every fluxRatio is relative to).
• models::Dict{String,model}: the per-line model object.
• lineCenters::Dict{String,Float64}: central wavelength per line, in any units – consistent units across lines in one CompositeModel are the caller's responsibility (only Δλ/λ matters downstream).
• fluxRatios::Dict{String,Float64}: each line's velocity-integrated line flux relative to lines[1] (so fluxRatios[lines[1]] == 1.0).

There is deliberately no lines(cm) accessor – cm.lines is the public API for the ordered name list (a lines(...) function would also collide with Workstream 5's lines::Tuple keyword).

CompositeModel(m::model; line::String, lineCenter::Float64)

Wrap an existing model as the first (reference) line of a new CompositeModel. The reference line's fluxRatio is fixed at 1.0 – every subsequent line's fluxRatio (see addLine!) is defined relative to it.

Parameters

• m::model: the model for this line (disk-wind, cloud, or a +-combined/raytraced model).
• line::String: unique name for this line (e.g. "Hα").
• lineCenter::Float64: central wavelength (any units, consistent across lines you plan to add to this CompositeModel). Must be positive.

Returns

• CompositeModel with one registered line.

rebuild(params::NamedTuple; overrides...) -> model
rebuild(m::model; overrides...) -> model

Reconstruct a model from the construction provenance recorded in params (see the model params field), optionally replacing individual constructor arguments via overrides.

rebuild dispatches on params.constructor:

• a leaf record (:DiskWindModel, :DiskWindModelMean, :cloudModel, :cloudModelVectors) calls the corresponding public constructor with its recorded arguments;
• a :+ record rebuilds left and right and sums them;
• a :raytrace! record rebuilds parent and re-raytraces it with the recorded raytrace arguments.

Overrides

overrides... broadcast to every leaf in which the key exists (so an override reaches all submodels of a combined/raytraced model that share that parameter). An override key that matches no leaf is an error (typo guard). Override semantics are exactly constructor semantics – whether an override moves geometry depends on the recorded parameterization. For example overriding α on a :DiskWindModelMean record moves the whole grid (α sets rMin/rMax via get_rMinMaxDiskWind), while overriding α on a :DiskWindModel (explicit rMin/rMax) record only reaches the intensity function – geometry is untouched. rebuild does exactly what calling the constructor with the merged arguments would do, nothing more.

Node-level arguments (IRatios, τCutOff, raytraceFreeClouds) are not override targets in v1: to change them, rebuild the parent and re-raytrace manually.

Reproducibility

Disk-wind construction is deterministic, so a rebuilt disk-wind line has identical geometry. Cloud models reproduce identical clouds only under rng=:philox with the recorded seed; a rng=:legacy record consumed the GLOBAL_RNG state, so its rebuild re-draws statistically-equivalent (different) gas from GLOBAL_RNG.

Not recorded

params captures construction history only. In-place mutators (removeNaN!, zeroDiskObscuredClouds!) are not recorded; a caller who ran them must re-apply the same calls in the same order to the rebuilt model. removeNaN! is observable-neutral (binned observables skip non-finite points), so the only unrecorded mutator with a physical effect is zeroDiskObscuredClouds!, which is deterministic and cheap to re-run.

rebuild(m::model; ...) is a convenience method forwarding to rebuild(m.params; ...); it errors when m.params === nothing (the model was not built by a public constructor).

DiskWindIntensity(;r::Union{Vector{Float64},Float64}, i::Float64, ϕ::Union{Vector{Float64},Float64},
        f1::Float64, f2::Float64, f3::Float64, f4::Float64, α::Float64, rMin::Float64 = 0.0, rMax::Float64 = Inf, _...)

Calculates the intensity from a disk-wind model of the BLR following the prescription given in Long+ 2023 and 2025, similar to Chiang and Murray 1996 and 1997, assumes optically thick line emission limit (Sobolev).

Arguments:

• r::Union{Vector{Float64},Float64} - radius from central mass (in terms of rₛ)
• i::Float64 - inclination angle (rad)
• ϕ::Union{Vector{Float64},Float64} - list of azimuthal angles (rad)
• f1::Float64 - strength of radial velocity gradient $\mathrm{d}v_{r}/\mathrm{d}r$
• f2::Float64 - strength of Keplerian shear velocity gradient $\mathrm{d}v_{\phi}/\mathrm{d}r$
• f3::Float64 - strength of radial lifting velocity gradient $\mathrm{d}v_{\theta}/\mathrm{d}r$
• f4::Float64 - strength of vertical lifting velocity gradient $\mathrm{d}v_{\theta}/\mathrm{d}\theta$ (equivalent to isotropic emission)
• α::Float64 - power-law index of the source function $S(r) \propto r^{-\alpha}$
• rMin::Float64 - minimum radius to consider (default: 0.0)
• rMax::Float64 - maximum radius to consider (default: Inf)

Returns:

• intensity (arbitrary units) as a Vector{Float64}

DiskWind_I_(r::Float64, ϕ::Float64, i::Float64, f1::Float64, f2::Float64, f3::Float64, f4::Float64, α::Float64)

Same as DiskWind_I_(r::Vector{Float64}, ϕ::Vector{Float64}, i::Float64, f1::Float64, f2::Float64, f3::Float64, f4::Float64, α::Float64) but for a single radius r and a single azimuthal angle ϕ.

DiskWind_I_(r::Float64, ϕ::Vector{Float64}, i::Float64, f1::Float64, f2::Float64, f3::Float64, f4::Float64, α::Float64)

Same as DiskWind_I_(r::Vector{Float64}, ϕ::Vector{Float64}, i::Float64, f1::Float64, f2::Float64, f3::Float64, f4::Float64, α::Float64) but for a single radius r and a vector of azimuthal angles ϕ.

DiskWind_I_(r::Vector{Float64}, ϕ::Float64, i::Float64, f1::Float64, f2::Float64, f3::Float64, f4::Float64, α::Float64)

Same as DiskWind_I_(r::Vector{Float64}, ϕ::Vector{Float64}, i::Float64, f1::Float64, f2::Float64, f3::Float64, f4::Float64, α::Float64) but for a vector of radii r and a single azimuthal angle ϕ.

DiskWind_I_(r::Vector{Float64}, ϕ::Vector{Float64}, i::Float64, f1::Float64, f2::Float64, f3::Float64, f4::Float64, α::Float64)

Calculates the intensity of the disk wind at radius r from the central mass and inclined at angle i (rad) over grid of azimuthal angles ϕ (rad) Follows the prescription given in Long+ 2023 and 2025, similar to Chiang and Murray 1996 and 1997, assumes optically thick line emission limit

Arguments:

• r::Vector{Float64} - radius from central mass (in terms of rₛ)
• ϕ::Vector{Float64} - list of azimuthal angles (rad)
• i::Float64 - inclination angle (rad)
• f1::Float64 - strength of radial velocity gradient dvᵣ/dr
• f2::Float64 - strength of Keplerian shear dvϕ/dr
• f3::Float64 - strength of velocity gradient dvθ/dr
• f4::Float64 - strength of velocity gradient dvθ/dθ (or isotropic emission)
• α::Float64 - power-law index of the source function S(r) ∝ r^(-α)

Returns:

• intensity (arbitrary units) as a Vector{Float64}

DiskWind_I_grid(r::Vector{Float64}, ϕ::Vector{Float64}, i::Float64, f1::Float64, f2::Float64,
    f3::Float64, f4::Float64, α::Float64, rMin::Float64, rMax::Float64)

Typed helper (function barrier) for DiskWindIntensity's vector path. Fills a single preallocated vector with the scalar DiskWind_I_ per point, then applies the (rMin, rMax) mask with the bounds rounded to 9 significant digits once instead of per point.

The outer mask guard replicates the original semantics exactly: when r contains NaN (out-of-grid sentinels), minimum(r)/maximum(r) are NaN, the guard comparisons are false, and no masking pass runs — NaN points keep their NaN intensity rather than being zeroed.

IsotropicIntensity(;r::Union{Vector{Float64},Float64}, ϕ::Union{Vector{Float64},Float64}, 
        rescale::Float64=1.0, rMin::Float64 = 0.0, rMax::Float64 = Inf, _...)

Returns a constant intensity value at all radii and azimuthal angles, rescaled by rescale factor.

IϕCloudMask(;r::Float64, ϕ::Float64, θₒ::Float64, ϕ₀::Float64, rot::Float64, i::Float64, κ::Float64=0.0, 
        ϕMin::Float64, ϕMax::Float64, overdense::Bool=false, _...)

Calculates the intensity using cloudIntensity but with an extra layer of masking based on azimuthal angle ϕ and the specified range [ϕMin, ϕMax]. If ϕ is outside this range, the intensity is set to 0.0. If overdense is true, the intensity is multiplied by 2.0 to account for the "lost" cloud and thus preserve total intensity in the model.

IϕDiskWindMask(;r::Union{Vector{Float64},Float64}, ϕ::Union{Vector{Float64},Float64}, i::Float64, 
        f1::Float64, f2::Float64, f3::Float64, f4::Float64, α::Float64, ϕMin::Float64, ϕMax::Float64, _...)

Calculates the intensity using DiskWindIntensity but with an extra layer of masking based on azimuthal angle ϕ and the specified range [ϕMin, ϕMax]. If ϕ is outside this range, the intensity is set to 0.0.

cloudIntensity(;r::Float64, ϕ::Float64, θₒ::Float64, ϕ₀::Float64, rot::Float64, i::Float64, κ::Float64=0.0, _...)

Calculate the intensity of the cloud at radius r from the central mass and inclined at angle i (rad) over grid of azimuthal angles ϕ (rad) following Pancoast+ 2011 and 2014.

Arguments

• r::Float64: radius from central mass (in terms of r₍ₛ₎)
• ϕ::Float64: azimuthal angle (rad)
• θₒ::Float64: opening angle of cloud
• ϕ₀::Float64: initial azimuthal angle
• rot::Float64: rotation angle
• i::Float64: inclination angle
• κ::Float64: anisotropy parameter

Returns

• Intensity (arbitrary units) as a Float64

vCirc(r::Float64, rₛ::Float64=1.0)

Calculate circular velocity at radius r from central mass, with Schwarzschild radius rₛ.

Defaults to rₛ=1.0 for unitless calculations.

vCircularCloud(;r::Float64, ϕ₀::Float64, i::Float64, rot::Float64, θₒ::Float64, rₛ::Float64=1.0, reflect::Bool=false, _...)

Calculate line of sight velocity for cloud in 3D space.

Arguments

• r::Float64: radius from central mass (in terms of rₛ)
• ϕ₀::Float64: starting azimuthal angle in ring plane (rad)
• i::Float64: inclination angle of ring plane (rad)
• rot::Float64: rotation of system plane about z axis (rad)
• θₒ::Float64: opening angle of point
• rₛ::Float64=1.0: Schwarzschild radius (optional, to convert to physical units)
• reflect::Bool=false: whether the point is reflected across the midplane of the disk
• _...: extra kwargs, ignored

Returns

• Line of sight velocity (Float64)

vCircularDisk(;r::Union{Float64,Vector{Float64}}, i::Float64, ϕ::Union{Vector{Float64},Float64}, rₛ=1.0, _...)

Calculate line of sight velocity for circular orbit at radius r from central mass and inclined at angle i (rad) over grid of azimuthal angles ϕ (rad).

vCircularRadialDisk(;r::Union{Float64,Vector{Float64}}, i::Float64, ϕ::Union{Vector{Float64},Float64}, 
        vᵣFrac::Union{Vector{Float64},Float64}=0.0, inflow::Union{Vector{Bool},Bool}=true, rₛ=1.0, _...)

Calculate line of sight velocity for circular orbit at radius r from central mass and inclined at angle i (rad) over grid of azimuthal angles ϕ (rad) with radial inflow/outflow.

vCloudTurbulentEllipticalFlow(;σρᵣ::Float64, σρc::Float64, σΘᵣ::Float64, σΘc::Float64, 
        θₑ::Float64, fEllipse::Float64, fFlow::Float64, σₜ::Float64, r::Float64, 
        i::Float64, rot::Float64, θₒ::Float64, rₛ::Float64=1.0, ϕ₀::Float64=0.0, 
        reflect::Bool=false, rng::AbstractRNG=Random.GLOBAL_RNG, _...)

Calculate line of sight velocity for cloud in 3D space with potential for elliptical orbital velocities, in/outflow, and turbulence as in Pancoast+14.

Arguments

• σρᵣ::Float64: Radial standard deviation around radial orbits
• σρc::Float64: Radial standard deviation around circular orbits
• σΘᵣ::Float64: Angular standard deviation around radial orbits
• σΘc::Float64: Angular standard deviation around circular orbits
• θₑ::Float64: Angle in vϕ-vr plane
• fEllipse::Float64: Fraction of elliptical orbits
• fFlow::Float64: If < 0.5, inflow, otherwise, outflow
• σₜ::Float64: Standard deviation of turbulent velocity
• r::Float64: Radius from central mass (in terms of rₛ)
• i::Float64: Inclination angle of ring plane (rad)
• rot::Float64: Rotation of system plane about z axis (rad)
• θₒ::Float64: Opening angle of point
• rₛ::Float64=1.0: Schwarzschild radius (optional, to convert to physical units)
• ϕ₀::Float64=0.0: Starting azimuthal angle in ring plane (rad)
• reflect::Bool=false: Whether the point is reflected across the midplane of the disk
• rng::AbstractRNG=Random.GLOBAL_RNG: Random number generator
• _...: Extra kwargs, ignored

Returns

• Line of sight velocity (Float64)

binAccumulate!(result, x, y, binEdges, overflow::Bool, uniform::Bool)

Typed accumulation loop for binnedSum (function barrier). When uniform is true the bin index is computed directly from the spacing (floor) and then corrected against the actual edges, which keeps the assignment bit-identical to searchsortedfirst — including the exact-edge convention (a point exactly on an interior edge belongs to the bin to its left) and any floating-point rounding of the direct computation. Non-uniform (user-supplied) edges use searchsortedfirst as before.

binModel(bins::Union{Int,Vector{Float64}}=100; m::model, yVariable::Union{String,Symbol,Function}, 
        xVariable::Union{String,Symbol,Function}=:v, kwargs...)
binModel(bins::Vector{Float64}, dx::Array{Float64,}; m::model, yVariable::Union{String,Symbol,Function}, 
        xVariable::Union{String,Symbol,Function}=:v, kwargs...)

Bin the model into a histogram, where each bin is the integrated value of the yVariable as a function of the xVariable.

Arguments

• m::model: Model object to bin
• yVariable::Union{String,Symbol,Function}: Variable to bin
  • Must be valid attribute of model.rings (e.g. :I, :v, :r, :i, :ϕ) or a function that can be applied to model.rings
  • Example: Keplerian disk time delays could be calculated like t(ring) = ring.r .* (1 .+ sin.(ring.ϕ) .* ring.i)
• bins::Union{Int,Vector{Float64}}: Number of bins or bin edges for binning
  • If Int: Number of bins with edges equally spaced between min/max of xVariable
  • If Vector{Float64}: Specific bin edges, with number of bins = length(bins)-1
  • Bins are left-edge exclusive, right-edge inclusive: a point exactly on an interior edge goes to the bin below it. Points at or below the first edge, or at or above the last edge, are only counted when overflow=true (into the first/last bin)
• xVariable::Union{String,Symbol,Function}=:v: Variable to bin over
  • Must be a valid attribute of model.rings or a function that can be applied to model.rings
• dx::Array{Float64,}: Integration element for each bin
  • If provided, used as the associated integral element
  • Otherwise defaults to ΔA in each ring struct

Returns

• Tuple{Vector{Float64},Vector{Float64},Vector{Float64}}: A tuple containing:
  • binEdges: Bin edges for the xVariable of the histogram
  • binCenters: Bin centers for the xVariable
  • yBinned: Binned values of the yVariables

binWeightedMean(m::model, x, yNum, yDen, bins, dx; kwargs...)

Shared helper for the weighted-mean profile options of getProfile (:delay, :r, :ϕ): bins yNum.*dA and yDen.*dA over x and returns (edges, centers, num./den). Replicates binModel's ΔA shape fix-up (per-ring scalar ΔA broadcast across a matrix y). Operates on precomputed arrays so the gathers can come from the model's memoized cache.

binnedSum(x::Array{Float64,}, y::Array{Float64, }; bins=100,
        overflow=false, centered=true, minX=nothing, maxX=nothing)

Bin the x and y variables into a histogram, where each bin is the sum of the y values for the corresponding x values.

Arguments

• x::Array{Float64,}: x variable to bin over
• y::Array{Float64,}: y variable to bin
• bins::Union{Int,Vector{Float64}}=100: Number of bins or bin edges for binning
  • If Int: Number of bins with edges equally spaced between min/max of x
  • If Vector{Float64}: Specific bin edges, with number of bins = length(bins)-1
  • Bins are left-edge exclusive, right-edge inclusive: a point exactly on an interior edge goes to the bin below it. Points at or below the first edge, or at or above the last edge, are only counted when overflow=true (into the first/last bin)
• overflow::Bool=false: If true, include values outside bin range in the first/last bins
• centered::Bool=true: If true, shift bin edges to center around middle value
• minX::Union{Float64,Nothing}=nothing: Minimum value of x for binning (defaults to minimum(x))
• maxX::Union{Float64,Nothing}=nothing: Maximum value of x for binning (defaults to maximum(x))

Returns

• Tuple{Vector{Float64},Vector{Float64},Vector{Float64}}: A tuple containing:
  • binEdges: Bin edges for the x variable
  • binCenters: Bin centers for the x variable
  • result: Binned sums of the y variable

binnedVariance(x::Array{Float64,}, θ::Array{Float64,}, w::Array{Float64,};
        bins=100, overflow=false, kwargs...)

Per-bin weighted variance of θ: bin the x variable as in binnedSum, then in each bin compute the w-weighted variance $\sigma^2 = \sum_k w_k(\theta_k-\mu)^2 / \sum_k w_k$ with μ the w-weighted mean of θ in that bin.

The variance is computed with a two-pass algorithm (pass 1: per-bin Σw and Σwθ → per-bin mean μ; pass 2: per-bin Σw(θ-μ)²/Σw) rather than the raw-moment identity σ² = ⟨θ²⟩ - ⟨θ⟩², which suffers catastrophic cancellation when the per-bin mean is large compared to the spread (e.g. a few clouds clustered tightly far from the bin mean). The two-pass form is also guaranteed non-negative. Bin membership in pass 2 follows the same rules as binnedSum, so the two can never disagree.

Arguments

• x::Array{Float64,}: x variable to bin over
• θ::Array{Float64,}: variable to take the variance of
• w::Array{Float64,}: weights (e.g. intensity × area)
• bins::Union{Int,Vector{Float64}}=100: number of bins or bin edges, as in binnedSum
• overflow::Bool=false: if true, include values outside the bin range in the first/last bins
• Additional kwargs (centered, minX, maxX) are passed to binnedSum

Returns

• Tuple{Vector{Float64},Vector{Float64},Vector{Float64}}: A tuple containing:
  • binEdges: Bin edges for the x variable
  • binCenters: Bin centers for the x variable
  • σ²: per-bin weighted variance of θ (NaN for empty bins)

constructBinEdges(bins::Int, minX::Float64, maxX::Float64, centered::Bool) -> Vector{Float64}

Build the uniform velocity-bin edges for an integer bin count, given the (NaN-resolved) data range minX/maxX. Shared by binnedSum and the GPU-resident profile path so the two cannot diverge on the centered edge convention.

getProfile(m::model, name::Union{String,Symbol,Function};
           bins::Union{Int,Vector{Float64}}=100,
           dx::Union{Array{Float64,},Nothing}=nothing, kwargs...)

Return a profile for the model based on the specified name.

Arguments

• m: Model object to get the profile from
• name: Name of the profile to get. Options include:
  • :line: Returns the line profile (integrated intensity as function of velocity)
  • :delay: Returns the delay profile (mean delays weighted by intensity/responsivity vs. velocity)
  • :r: Returns the mean radius (weighted by intensity) as function of velocity
  • :ϕ: Returns the mean azimuthal angle (weighted by intensity) as function of velocity
  • :phase: Returns the phase profile (integrated phase as function of velocity)
    • Requires U [Mλ], V [Mλ], PA [rad], and BLRAng [rad] as keyword arguments
  • :moment2: Returns the per-channel squared angular size of the image along the baseline direction(s) [rad²]
    • Requires U [Mλ], V [Mλ], PA [rad], and BLRAng [rad] as keyword arguments
    • Use secondMoment directly to also obtain the line-integrated size (a scalar per baseline, which does not fit in a profile struct)
  • Function: Returns the intensity weighted mean of this function vs. velocity
• bins: Number of bins or bin edges for binning
  • If Int: Number of bins with edges equally spaced between min/max velocity
  • If Vector{Float64}: Specific bin edges, with number of bins = length(bins)-1
  • Bins are left-edge exclusive, right-edge inclusive: a point exactly on an interior edge goes to the bin below it. Points at or below the first edge, or at or above the last edge, are only counted when overflow=true (into the first/last bin)
• dx: Integration element for each ring (defaults to ΔA in each ring struct if nothing)
• Additional kwargs passed to binModel include:
  • overflow=true: Include overflow bins
  • centered=true: Center bins around 0.0
  • minX, maxX: Set min/max bin boundaries

Returns

• profile: A profile object containing bin edges, bin centers, and binned sums. Assign to a model object using setProfile! to store the profile in the model.

A device-resident method (getProfile(::ResidentModel, …)) runs this on the GPU without re-flattening the model each call; see gpu and ResidentModel.

phase(m::model; U, V, PA, BLRAng, returnAvg=false, offAxisInds=nothing, kwargs...)

Calculate differential phase signal for a model based on specified baselines, model orientation, and BLR angular size.

Arguments

• m::model: Model object to calculate phase for
• U::Vector{Float64}: U component of complex visibility in [Mλ]
• V::Vector{Float64}: V component of complex visibility in [Mλ]
• PA::Float64: On-sky position angle of the model in radians
• BLRAng::Float64: Characteristic size of the BLR model in radians (conversion from $r_s$ to radians)
• returnAvg::Bool=false: If true, returns the average phase across all baselines
• offAxisInds::Union{Nothing,Vector{Int}}=nothing: If provided, only calculates phase for baselines at specified indices

Returns

• If returnAvg=true: Tuple{Vector{Float64},Vector{Float64},Vector{Float64}} containing:
  • Bin edges for velocity
  • Bin centers for velocity
  • Average differential phase (in radians)
• If returnAvg=false: Vector{Tuple{Vector{Float64},Vector{Float64},Vector{Float64}}} containing:
  • For each baseline, a tuple of bin edges, bin centers, and differential phase

The differential phase is calculated by integrating the phase over the model at each velocity bin, weighted by the intensity and area of each ring element.

A device-resident method (phase(::ResidentModel, …)) runs this on the GPU without re-flattening the model each call; see gpu and ResidentModel.

secondMoment(m::model; U, V, PA, BLRAng, returnAvg=false, offAxisInds=nothing, kwargs...)

Calculate the second central moment (squared angular size) of the model image projected along interferometric baseline directions, both per velocity channel and integrated over the whole line.

For each baseline the per-channel size is the intensity-weighted variance of the projected position û⋅x in each velocity bin, $\sigma^2_\theta(v, \hat{u}) = \langle (\hat{u}\cdot(x - \bar{x}_v))^2 \rangle_v$, i.e. the squared rms angular extent of the channel image along the baseline direction. In the marginally-resolved limit this predicts the differential visibility-amplitude dip via $|V_l|(v) = 1 - 2\pi^2 |u|^2 \sigma^2_\theta(v)$ with $|u| = \sqrt{U^2+V^2}\times10^6$ (baseline length in units of the observed wavelength).

The line-integrated size is the same variance taken over the whole line image (all velocities in a single channel), which by the law of total variance equals the flux-weighted mean of the per-channel sizes plus the flux-weighted spread of the per-channel photocenters — for a rotation-dominated BLR the photocenter-spread term is comparable to the mean per-channel size, so the line-integrated size is significantly larger than the flux-averaged channel size.

Note that these are moments of the line image only: unlike phase, no continuum dilution factor f/(1+f) is applied.

Arguments

• m::model: Model object to calculate image second moments for
• U::Vector{Float64}: U component of complex visibility in [Mλ]
• V::Vector{Float64}: V component of complex visibility in [Mλ]
• PA::Float64: On-sky position angle of the model in radians
• BLRAng::Float64: Characteristic size of the BLR model in radians (conversion from $r_s$ to radians)
• returnAvg::Bool=false: If true, returns the average across all baselines
• offAxisInds::Union{Nothing,Vector{Int}}=nothing: If provided, only calculates moments for baselines at specified indices
• Additional kwargs (bins, overflow, centered, minX, maxX) control the velocity binning as in binnedSum

Returns

• If returnAvg=true: Tuple{Vector{Float64},Vector{Float64},Vector{Float64},Float64} containing:
  • Bin edges for velocity
  • Bin centers for velocity
  • Average per-channel squared angular size σ²(v) (in rad², NaN for empty bins)
  • Average line-integrated squared angular size σ²tot (in rad²)
• If returnAvg=false: Vector{Tuple{Vector{Float64},Vector{Float64},Vector{Float64},Float64}} containing:
  • For each baseline, a tuple of bin edges, bin centers, per-channel σ²(v), and line-integrated σ²tot

A device-resident method (secondMoment(::ResidentModel, …)) runs this on the GPU without re-flattening the model each call; see gpu and ResidentModel.

setProfile!(m::model, p::profile; overwrite::Bool=false)

Set a profile in the model's profiles dictionary. If the profile already exists and overwrite is false, a warning is issued.

t(ring::ring, subFxn::Function=tDisk)

Calculate time delays for a point in a disk or cloud using a custom function subFxn that takes a ring struct.

It is more peformant to pass the function directly rather than figure it out on the fly if known ahead of time.

t(ring::ring)

Calculate time delays for a point in a disk as `` t = \eta r \left(1 + \cos(\phi) \sin(i)\right)`` or a cloud with opening angle ``\theta_o``
as the x-coordinate of the point subtracted from the radial distance of the point ``t = r - x``.

tCloud(ring::ring)

Calculate time delays from the general geometric formula ``t = \eta (r - x)``, where ``x`` is the
3D system coordinate of the point towards the camera (cached on the ring, see `getXYZ`).

Despite the name this is exact for *any* geometry: for a flat ring (``\theta_o = 0``) it reduces
algebraically to the `tDisk` formula (``x = r\cos\phi\sin i``), and unlike `tDisk` it also
accounts for lifted (``\theta_o \neq 0``) and reflected points. Works on both point (scalar) and
continuous (vector) rings.

tDisk(ring::ring)

Calculate time delays for a point in a disk as `` t = \eta r \left(1 - \cos(\phi) \sin(i)\right)``

getΨ(m::model,vEdges::Array{Float64},tEdges::Array{Float64})

Calculate the 2D transfer function Ψ for a model m over specified velocity and time bins, whose edges are given by vEdges and tEdges.

A device-resident method (getΨ(::ResidentModel, …)) runs this on the GPU without re-flattening the model each call; see gpu and ResidentModel.

getΨ(m::model,vBins::Int64,tBins::Int64)

Calculate the 2D transfer function Ψ for a model m over specified number of velocity bins vBins and time bins tBins. The velocity and time edges are automatically calculated based on the minimum and maximum values for velocity and delays in the model.

getΨt(m::model,tEdges::Array{Float64},overflow::Bool=false;)

Calculate the 1D transfer function Ψ(t) for a model m over specified time edges tEdges. The overflow parameter determines whether to include contributions from delays outside the specified edges in the edge bins.

A device-resident method (getΨt(::ResidentModel, …)) runs this on the GPU without re-flattening the model each call; see gpu and ResidentModel.

getΨt(m::model,tBins::Int64,maxT::Float64=Inf,overflow::Bool=false)

Calculate the 1D transfer function Ψ(t) for a model m over specified number of time bins tBins. The maxT parameter specifies the maximum time delay to consider, and overflow determines whether to include contributions from delays outside the specified edges in the edge bins.

response(r::Float64, ηₒ::Float64, η₁::Float64, αRM::Float64, rNorm::Float64)

Positional-argument version of response for hot loops — avoids the per-call keyword splat.

response(r::Float64; ηₒ::Float64=0.5, η₁::Float64=0.5, αRM::Float64=0.0, rNorm::Float64=1.0, _...)

Calculate response function for use in reverberation mapping calculations.

ΨAccumulate!(Ψ, v, delays, I, ΔA, vEdges, tEdges)

Single-pass accumulation for getΨ (function barrier): one sweep over all points, locating each point's (velocity, delay) bin instead of building a boolean mask per bin pair (the old approach was O(nbins² × npoints) in time and memory traffic).

Bin assignment uses the same convention as binnedSum (searchsortedfirst): interior edges are left-EXCLUSIVE / right-inclusive — a point exactly on an interior edge goes to the bin below it — and points at or below the first edge or at or above the last edge are dropped (Ψ has no overflow bins). NaN velocities/delays never land in a bin; NaN intensities/areas still accumulate into their bin so the final floor maps the poisoned bin to 1e-30. This matches the line profile and every other binnedSum-based observable, so a shared edge set bins identically everywhere.

ΨtAccumulate!(Ψt, delays, I, ΔA, tEdges) -> (sUnder, sOver)

Single-pass accumulation for getΨt (function barrier – the gathered arrays are not inferrable in the caller). Same binnedSum (searchsortedfirst, left-exclusive interior edge) convention as ΨAccumulate!; additionally collects the underflow/overflow sums (delays <= tEdges[1] / delays >= tEdges[end], NaN in neither) in the same sweep, matching binnedSum's overflow=true edge-bin folding.

photograph(r::Float64, ϕ₀::Float64, i::Float64, rot::Float64, θₒ::Float64, reflect::Bool=false)

Calculate the image coordinates from system coordinates r, ϕ + inclination angle i.

Arguments

• r::Float64: radius from central mass (in terms of rₛ)
• ϕ₀::Float64: unrotated azimuthal angle in ring plane (rad)
• i::Float64: inclination angle of ring plane (rad)
• rot::Float64: rotation of system plane about z axis (rad)
• θₒ::Float64: ring opening angle
• reflect::Bool=false: whether the point is reflected across the midplane of the disk

Returns

• α::Float64: image x coordinate (in terms of rₛ)
• β::Float64: image y coordinate (in terms of rₛ)

Note

This function is coordinate photography only. To visualize models, see Image.`

raytrace!(rm::ResidentModel; IRatios=1.0, τCutOff=1.0, raytraceFreeClouds=false) -> ResidentModel

Device-resident dispatch of raytrace!: run the full bin→sort→segmented-scan→ compact combine on rm's backend with no host round-trips, returning a freshly raytraced ResidentModel (drop-in for the observable methods). rm must carry raytrace metadata — build it with gpu(m) (CUDA loaded) or resident(m; raytrace=true).

Like raytrace!(::model) this returns a new combined handle rather than mutating in place (the point count changes: one combined point per active output pixel + surviving free clouds), so write rm = raytrace!(rm). Semantics match the host raytrace!: IRatios is a global per-submodel weight, τCutOff the optical-depth cutoff, raytraceFreeClouds enables cloud–cloud attenuation.

raytrace!(m::model; IRatios::Union{Float64,Array{Float64,}}=1.0,
        τCutOff::Float64=1.0, raytraceFreeClouds::Bool=false,
        backend=nothing, T=Float64)

Perform raytracing for a model, combining overlapping components along line of sight.

This function should be called after combining all relevant models (i.e. mCombined = m1 + m2 + m3...). Points are binned by pixel, sorted front-to-back along the line of sight, and combined with a segmented exp(-τ) scan: each point attenuates those behind it, points past τCutOff are dropped, and IRatios is applied as a global per-submodel weight. A new model object is returned with the extraneous points removed. Note that this function will mutate the input model objects.

Arguments

• m::model: Model to raytrace
• IRatios::Union{Float64,Array{Float64,}}=1.0: Global emissivity weights for each submodel
  • If Float64, applies to all submodels equally
  • If array, applies to each submodel individually (must match number of submodels)
  • Used when combining models with different intensity functions if they aren't properly normalized
• τCutOff::Float64=1.0: Maximum optical depth to raytrace to (stops when τ > τCutOff)
• raytraceFreeClouds::Bool=false: Whether to raytrace free clouds (cloud-cloud raytracing)
  • If false, clouds are only raytraced if they overlap with a continuous model
  • If true, clouds will be checked for overlap with other clouds and raytraced accordingly
• backend=nothing: Optional KernelAbstractions backend to run the bin→sort→scan on a device.
  • If nothing, uses the CPU implementation.
  • Pass CUDABackend() (with CUDA.jl loaded) to run on the GPU; see also gpu.
• T=Float64: Element type for the device arrays when backend is set (use Float32 on GeForce GPUs for speed).

Returns

• m::model: Model with raytraced points

raytrace(α::Float64, β::Float64, i::Float64, rot::Float64, θₒPoint::Float64)

Calculate where ray traced back from camera coordinates α and β intersects the system (assumes circular geometry).

Arguments

• α::Float64: image x coordinate (in terms of rₛ)
• β::Float64: image y coordinate (in terms of rₛ)
• i::Float64: inclination angle of system (rad)
• rot::Float64: how the point was rotated about z axis (rad)
• θₒPoint::Float64: opening angle of current point

Returns

• r::Float64: distance from central mass (in terms of rₛ)
• ϕ::Float64: azimuthal angle of system ring plane at intersection
• ϕ₀::Float64: original azimuthal angle in ring plane (no rotation)

Note

This function is coordinate raytracing only. To raytrace models and combine intensities, see raytrace!.

raytrace(α::Float64, β::Float64, i::Float64, rot::Float64, θₒPoint::Float64, M::Matrix{Float64}, r3D::Matrix{Float64})

Same as the 6-argument raytrace(α, β, i, rot, θₒPoint, M) method but additionally returns the 3D system coordinates of the intersection point (camera at +x), computed from the supplied r3D = get_r3D(i, rot, θₒPoint) matrix. Used at model construction to fill the per-point coordinate cache (ring.x/y/z) without a second rotation pass.

Returns

• r, ϕ, ϕ₀, x, y, z (all Float64)

raytrace(α::Float64, β::Float64, i::Float64, rot::Float64, θₒPoint::Float64, 
        r3D::Matrix{Float64}, xyz::Vector{Float64}, matBuff::Matrix{Float64}, 
        colBuff::Vector{Float64})

Performant version of raytrace function – calculate where ray traced back from camera coordinates α, β intersects the system (assumes circular geometry).

Arguments

• α::Float64: image x coordinate (in terms of rₛ)
• β::Float64: image y coordinate (in terms of rₛ)
• i::Float64: inclination angle of system (rad)
• rot::Float64: rotation of current point about z axis (rad)
• θₒPoint::Float64: opening angle of current point
• r3D::Matrix{Float64}: matrix that rotates system plane into XY plane
• xyz::Vector{Float64}: preallocated xyz vector (but not precalculated)
• matBuff::Matrix{Float64}: preallocated buffer matrix for storing result of 3x3 matrix multiplication
• colBuff::Vector{Float64}: preallocated buffer vector for storing final matrix multiplication result

Returns

• r::Float64: distance from central mass (in terms of rₛ)
• ϕ::Float64: azimuthal angle of system ring plane at intersection
• ϕ₀::Float64: original azimuthal angle in ring plane

Note

This function is coordinate raytracing only. To raytrace models and combine intensities, see raytrace!.

raytrace(α::Float64, β::Float64, i::Float64, rot::Float64, θₒPoint::Float64, M::Matrix{Float64})

Allocation-free version of raytrace — calculate where ray traced back from camera coordinates α, β intersects the system (assumes circular geometry), using a single precomputed rotation matrix.

M must be undo_tilt * get_r3D(i, rot, θₒPoint) where undo_tilt = [sin(i) 0.0 -cos(i); 0.0 1.0 0.0; cos(i) 0.0 sin(i)] — both factors are constant for fixed (i, rot, θₒ) so M should be computed once outside any pixel loop.

Returns

• r::Float64: distance from central mass (in terms of rₛ)
• ϕ::Float64: azimuthal angle of system ring plane at intersection
• ϕ₀::Float64: original azimuthal angle in ring plane

Note

Supersedes the buffer-based 9-argument method (which is retained for compatibility).

removeDiskObscuredClouds!(m::model, rotate3D::Function=rotate3D_scalar)

Remove clouds that are obscured by the disk.

Performs simple raytracing for an optically thick obscuring disk. The function modifies the input model by removing cloud points that are obscured by the disk. Note that this is a mutating operation and the input model will be modified in place.

Arguments

• m::model: Model to remove disk obscured clouds. Should be a combined model consisting of a disk component and a cloud component.
• rotate3D::Function=rotate3D_scalar: Function to rotate coordinates in 3D space (must return an indexable (x,y,z))

Returns

• m::model: Model with disk obscured clouds removed

See also

• zeroDiskObscuredClouds!: Function to zero out disk obscured clouds instead of removing them

zeroDiskObscuredClouds!(m::model; diskCloudIntensityRatio::Float64=1.0, rotate3D::Function=rotate3D_scalar)

Zero out the intensities of clouds that are obscured by the disk.

Performs simple raytracing for an optically thick obscuring disk. The function modifies the input model by setting the intensity of obscured cloud points to zero and adjusting the disk intensity according to the specified ratio.

Arguments

• m::model: Model to zero out disk obscured clouds. Should be a combined model consisting of a disk component and a cloud component.
• diskCloudIntensityRatio::Float64=1.0: Ratio of disk to cloud intensity, used to scale the disk intensities after zeroing out clouds
• rotate3D::Function=rotate3D_scalar: Function to rotate coordinates in 3D space (must return an indexable (x,y,z))

Returns

• m::model: Model with disk obscured clouds zeroed out

See also

• removeDiskObscuredClouds!: Function to remove disk obscured clouds instead of zeroing them out

RaytraceMeta

Device-resident metadata that raytrace!(::ResidentModel) needs beyond the flat ModelArrays: a list of output-grid blocks (grids, one per continuous submodel — each a _RtGridBlock with grid edges, per-output-pixel area/camera coords, and the block's inner radius for ordering) plus the per-point submodel index and discrete/cloud flag (submodel/discrete, length = number of input points). raytrace! assembles the blocks innermost-first into a single output grid at call time, so any number/kind of submodels combine generically (mirroring the host _rt_build_output). All array fields live on the same backend as the model columns.

addGrid!(m::model, colors=nothing, nϕ::Int=64)

Add a grid to the model image plot - mostly a debugging tool to visualize grid cells of overlapping models.

Arguments

• m::model: Model object to add grid to
• colors=nothing: Vector of colors for each submodel (if nothing, uses default colors)
• nϕ::Int=64: Number of azimuthal angles to use for the grid

Returns

• ::Plots.plot: Plot with grid added

expandPerPoint(m::model, variable::Symbol)

Like getVariable(m, variable; flatten=true) but per-ring scalar values are expanded to per-point length, so the result always aligns elementwise with getVariable(m, :I, flatten=true) — including for combined models where some rings store a field (e.g. η) as a single scalar while I is a vector.

The within-chunk ordering matches getVariable's flattening (vec(stack(chunk, dims=1)), i.e. point-major across the rings of each submodel), so elementwise products with other flattened gathers are valid. Results are memoized in m.cache under (variable, :perpoint).

gatherStack(rings::Vector{ring}, variable::Symbol)

Typed helper (function barrier) for getVariable: stacks getfield(ring, variable) over all rings into a Vector (one scalar per ring, e.g. cloud models) or a Matrix with one row per ring (vector-valued fields), without the intermediate vector-of-vectors that stack requires. Preserves the element type of the field (e.g. Bool for :reflect).

getFlattenedCameraIndices(m::model)

Get flattened camera indices corresponding to rings in model.

Arguments

• m::model: Model object to extract camera indices from

Returns

• camStartInds::Vector{Int}: Vector of camera starting indices with length equal to m.subModelStartInds

getRingFromFlattenedInd(m::model, flattenedInd::Int) -> Tuple{Int, Int}

Retrieve the model ring index and subindex (if applicable) from flattened array index.

Arguments

• m::model: Model object with rings
• flattenedInd::Int: The index in the flattened array we need to work back from

Returns

• row::Int: The ring index in model.rings that the flattened index corresponds to
• column::Int: The subindex that matches the flattened index passed to this function.

getVariable(m::model, variable::Function; flatten=false)

Retrieve model variable when specified as a Function. See main docstring for details.

Note

Results are memoized in m.cache only for the package's own pure delay functions (t, tDisk, tCloud) — user-supplied functions are always re-evaluated, since they may capture state the cache cannot see.

getVariable(m::model, variable::Union{String,Symbol,Function}; flatten=false)

Retrieve elements from model object and stack them into matrices for easy manipulation.

Arguments

• m::model: Model object to extract variables from
• variable::Union{String,Symbol,Function}: Variable to extract from model
  • If String, will be converted to Symbol
  • Must be a valid attribute of model.rings (e.g. :I, :v, :r, :i, :ϕ) or a function that can be applied to model.rings
  • Example: Keplerian disk time delays could be calculated like t(ring) = ring.r .* (1 .+ sin.(ring.ϕ) .* ring.i)
• flatten::Bool=false: If true, flatten the result to a vector

Returns

• Array{Float64,}: Matrix or vector of extracted variable from model.rings, created by stacking the output variable for each ring
  • For example, if variable given is :I, result will have shape (length(r), length(ϕ)) as at each r and ϕ there is a value of I
  • If flatten=true, result will be a flattened vector

getVariable(m::model, variable::Symbol; flatten=false)

Retrieve model variable when specified as a Symbol. See main docstring for details.

getXYZ(ring::ring)

Return the 3D system coordinates (x, y, z) of every point in ring (camera at +x).

Uses the cached ring.x/y/z fields when present; otherwise computes them with rotate3D_scalar and stores them on the (mutable) ring so the transform is only ever done once per point.

Returns

• For a cloud/point ring (scalar fields): NTuple{3,Float64}
• For a continuous ring (vector fields): Tuple{Vector{Float64},Vector{Float64},Vector{Float64}}

get_r3D(i::Float64, rot::Float64, θₒ::Float64) -> Matrix{Float64}

Calculate rotation matrix to transform from initial XY plane coordinates to 3D space.

Parameters

• i::Float64: Inclination angle of ring (rad)
• rot::Float64: Rotation of ring plane about z axis (rad)
• θₒ::Float64: Opening angle of point (rad)

Returns

• matrix::Matrix{Float64}: 3×3 rotation matrix

get_rMinMaxDiskWind(r̄::Float64, rFac::Float64, α::Float64)

Calculate the minimum and maximum radius of model given the (intensity weighted) mean radius r̄, the radius factor rFac, and the power-law index α following Long+ 2023.

Parameters

• r̄::Float64: mean radius of model (in terms of rₛ)
• rFac::Float64: radius scaling factor
• α::Float64: power-law index of the source function $S(r) \propto r^{-\alpha}$ (cannot be 1/2 or 3/2 as this divides by zero)

Returns

• rMin::Float64: minimum radius of model (in terms of rₛ)
• rMax::Float64: maximum radius of model (in terms of rₛ)

BLR.image(m::model, variable::Union{String,Symbol,Function}, kwargs...)

Generate an image of the model where the color of each point is determined by the variable provided.

Arguments

• m::model: Model object to extract variable from
• variable::Union{String,Symbol,Function}: Variable to extract from model
  • If String, will be converted to Symbol
  • Must be a valid attribute of model.rings (e.g. :I, :v, :r, :i, :ϕ) or a function that can be applied to model.rings
  • Example: Keplerian disk time delays could be calculated like t(ring) = ring.r .* (1 .+ sin.(ring.ϕ) .* ring.i)

Keywords

• Additional keyword arguments are passed to Plots.plot

Returns

• p::Plots.plot: Plot object representing the generated image

BLR.image(m::model, variable::Union{String,Symbol,Function}, kwargs...)

Generate an image of the model where the color of each point is determined by the variable provided.

Arguments

• m::model: Model object to extract variable from
• variable::Union{String,Symbol,Function}: Variable to extract from model
  • If String, will be converted to Symbol
  • Must be a valid attribute of model.rings (e.g. :I, :v, :r, :i, :ϕ) or a function that can be applied to model.rings
  • Example: Keplerian disk time delays could be calculated like t(ring) = ring.r .* (1 .+ sin.(ring.ϕ) .* ring.i)

Keywords

• Additional keyword arguments are passed to Plots.plot

Returns

• p::Plots.plot: Plot object representing the generated image

BLR.image(m::model, variable::Union{String,Symbol,Function}, kwargs...)

Generate an image of the model where the color of each point is determined by the variable provided.

Arguments

• m::model: Model object to extract variable from
• variable::Union{String,Symbol,Function}: Variable to extract from model
  • If String, will be converted to Symbol
  • Must be a valid attribute of model.rings (e.g. :I, :v, :r, :i, :ϕ) or a function that can be applied to model.rings
  • Example: Keplerian disk time delays could be calculated like t(ring) = ring.r .* (1 .+ sin.(ring.ϕ) .* ring.i)

Keywords

• Additional keyword arguments are passed to Plots.plot

Returns

• p::Plots.plot: Plot object representing the generated image

midPlaneXZ(x::Float64, i::Float64) -> Float64

Height of the disk mid-plane in the x–z plane at camera coordinate x for inclination i; used to decide cloud reflection.

plot3d(m::model, [variable], [annotate], [kwargs...])
plot3d(rm::ResidentModel, [variable], [annotate], [kwargs...])

Generate a 3D plot of the model geometry, optionally colored by a variable. Accepts a host model or a ResidentModel (its flat columns are copied to the host; the cached x/y/z are used directly).

Parameters

• m::model / rm::ResidentModel: Model object to plot
• variable::Union{String,Symbol,Function}=nothing: Variable to color the points by
  • If String, will be converted to Symbol
  • Must be a valid attribute of model.rings (e.g., :I, :v, :r, :i, :ϕ) or a function that can be applied to model.rings
  • Example: Keplerian disk time delays could be calculated like t(ring) = ring.r .* (1 .+ sin.(ring.ϕ) .* ring.i)
  • If not provided, defaults to nothing (no coloring). For a ResidentModel, a Function variable is not supported (pass a column Symbol); a combined (multi-submodel) ResidentModel must carry raytrace metadata (build it with raytrace=true) so the series can be split per submodel.
• annotate::Bool=true: Whether to annotate the camera position and model orientation in the plot
• kwargs...: Additional keyword arguments passed to Plots.plot

Returns

• p::Plots.plot: 3D plot of the model geometry, optionally colored by the variable provided

plot3d(m::model, [variable], [annotate], [kwargs...])
plot3d(rm::ResidentModel, [variable], [annotate], [kwargs...])

Generate a 3D plot of the model geometry, optionally colored by a variable. Accepts a host model or a ResidentModel (its flat columns are copied to the host; the cached x/y/z are used directly).

Parameters

• m::model / rm::ResidentModel: Model object to plot
• variable::Union{String,Symbol,Function}=nothing: Variable to color the points by
  • If String, will be converted to Symbol
  • Must be a valid attribute of model.rings (e.g., :I, :v, :r, :i, :ϕ) or a function that can be applied to model.rings
  • Example: Keplerian disk time delays could be calculated like t(ring) = ring.r .* (1 .+ sin.(ring.ϕ) .* ring.i)
  • If not provided, defaults to nothing (no coloring). For a ResidentModel, a Function variable is not supported (pass a column Symbol); a combined (multi-submodel) ResidentModel must carry raytrace metadata (build it with raytrace=true) so the series can be split per submodel.
• annotate::Bool=true: Whether to annotate the camera position and model orientation in the plot
• kwargs...: Additional keyword arguments passed to Plots.plot

Returns

• p::Plots.plot: 3D plot of the model geometry, optionally colored by the variable provided

plot3d(m::model, [variable], [annotate], [kwargs...])
plot3d(rm::ResidentModel, [variable], [annotate], [kwargs...])

Generate a 3D plot of the model geometry, optionally colored by a variable. Accepts a host model or a ResidentModel (its flat columns are copied to the host; the cached x/y/z are used directly).

Parameters

• m::model / rm::ResidentModel: Model object to plot
• variable::Union{String,Symbol,Function}=nothing: Variable to color the points by
  • If String, will be converted to Symbol
  • Must be a valid attribute of model.rings (e.g., :I, :v, :r, :i, :ϕ) or a function that can be applied to model.rings
  • Example: Keplerian disk time delays could be calculated like t(ring) = ring.r .* (1 .+ sin.(ring.ϕ) .* ring.i)
  • If not provided, defaults to nothing (no coloring). For a ResidentModel, a Function variable is not supported (pass a column Symbol); a combined (multi-submodel) ResidentModel must carry raytrace metadata (build it with raytrace=true) so the series can be split per submodel.
• annotate::Bool=true: Whether to annotate the camera position and model orientation in the plot
• kwargs...: Additional keyword arguments passed to Plots.plot

Returns

• p::Plots.plot: 3D plot of the model geometry, optionally colored by the variable provided

profile(m::model, [variable], [kwargs...])
profile(rm::ResidentModel, variable, [kwargs...])
profile(cm::CompositeModel, [name], [kwargs...])
profile(rcm::ResidentCompositeModel, [name], [kwargs...])
profile(p::profile, [kwargs...])

Plot line profiles, normalized to the maximum value of each profile.

Arguments

• model: A model object containing profile data, or a ResidentModel.
• variable: A symbol or string (or list of symbols/strings) specifying which profile to plot. For a host model this is optional — if omitted, every profile already set on the model is plotted. A ResidentModel stores no preset profiles, so a variable is required there and the profile is computed on demand via getProfile (choose from :line, :delay, :r, :ϕ, :phase, :moment2).
• kwargs...: Additional keyword arguments passed to Plots.plot

CompositeModel and bare profile support (Workstream 4/5)

• cm::CompositeModel: plots the requested name profile (default :line, computed on demand via getProfile(cm, name; line=...)) for every line, one series per line labeled by the line name and normalized by that line's own max |binSums| (the same convention the multi-profile model case below uses) – so lines are shape-comparable regardless of fluxRatio. The one exception is profile(cm, :ratio, (a, b)) (W5): the velocity-resolved line ratio is computed from the PAIR of lines via getProfile(cm, :ratio; lines=(a, b)) and rendered as a single unnormalized series (the ratio's absolute scale is the point), identical to plotting the returned bare profile. The pair is a positional argument because Plots reserves lines as a magic line-attribute alias – a lines=... keyword would be consumed by the attribute preprocessing and never reach the recipe. Other :ratio keywords (bins, floor, ...) are not forwarded – for those, call getProfile yourself and plot the returned struct.
• rcm::ResidentCompositeModel: same as the CompositeModel case (one normalized series per line), except each line's profile is computed on its device-resident columns via getProfile(rcm, name; line=...) (W4-G3), with the resident restrictions (uniform bins, no custom Function profiles).
• p::profile: plots a single bare profile struct directly (binCenters vs binSums, unnormalized). Needed so Workstream 5's :ratio profile (also a plain profile struct) is directly plottable. Note profile is both this struct's type and (via @userplot Profile above) the name of the generated plotting function – they're the same generic function, but the @kwdef constructor takes zero positional arguments (kwargs only), while this plotting method takes one-or-more positional arguments, so profile(p) is never ambiguous with profile(name=..., binCenters=..., ...).

profile(m::model, [variable], [kwargs...])
profile(rm::ResidentModel, variable, [kwargs...])
profile(cm::CompositeModel, [name], [kwargs...])
profile(rcm::ResidentCompositeModel, [name], [kwargs...])
profile(p::profile, [kwargs...])

Plot line profiles, normalized to the maximum value of each profile.

Arguments

• model: A model object containing profile data, or a ResidentModel.
• variable: A symbol or string (or list of symbols/strings) specifying which profile to plot. For a host model this is optional — if omitted, every profile already set on the model is plotted. A ResidentModel stores no preset profiles, so a variable is required there and the profile is computed on demand via getProfile (choose from :line, :delay, :r, :ϕ, :phase, :moment2).
• kwargs...: Additional keyword arguments passed to Plots.plot

CompositeModel and bare profile support (Workstream 4/5)

• cm::CompositeModel: plots the requested name profile (default :line, computed on demand via getProfile(cm, name; line=...)) for every line, one series per line labeled by the line name and normalized by that line's own max |binSums| (the same convention the multi-profile model case below uses) – so lines are shape-comparable regardless of fluxRatio. The one exception is profile(cm, :ratio, (a, b)) (W5): the velocity-resolved line ratio is computed from the PAIR of lines via getProfile(cm, :ratio; lines=(a, b)) and rendered as a single unnormalized series (the ratio's absolute scale is the point), identical to plotting the returned bare profile. The pair is a positional argument because Plots reserves lines as a magic line-attribute alias – a lines=... keyword would be consumed by the attribute preprocessing and never reach the recipe. Other :ratio keywords (bins, floor, ...) are not forwarded – for those, call getProfile yourself and plot the returned struct.
• rcm::ResidentCompositeModel: same as the CompositeModel case (one normalized series per line), except each line's profile is computed on its device-resident columns via getProfile(rcm, name; line=...) (W4-G3), with the resident restrictions (uniform bins, no custom Function profiles).
• p::profile: plots a single bare profile struct directly (binCenters vs binSums, unnormalized). Needed so Workstream 5's :ratio profile (also a plain profile struct) is directly plottable. Note profile is both this struct's type and (via @userplot Profile above) the name of the generated plotting function – they're the same generic function, but the @kwdef constructor takes zero positional arguments (kwargs only), while this plotting method takes one-or-more positional arguments, so profile(p) is never ambiguous with profile(name=..., binCenters=..., ...).

reflect!(xyzSys::Vector{Float64}, i::Float64) -> Vector{Float64}

Reflect coordinates in 3D space across the ring plane.

Parameters

• xyzSys::Vector{Float64}: [x;y;z] coordinates in 3D space
• i::Float64}: Inclination angle of ring plane (rad)

Returns

• xyzSys::Vector{Float64}: [x';y';z'] coordinates in 3D space after reflection

reflect_scalar(x::Float64, y::Float64, z::Float64, i::Float64) -> NTuple{3,Float64}

Allocation-free scalar version of reflect! — reflect coordinates in 3D space across the ring plane. Same math as reflect! but takes/returns plain scalars instead of mutating a Vector.

removeNaN!(m::model)

Remove points with I = NaN from model.

Parameters

• m::model: Model to remove points from

Returns

• m::model: Model with NaN points removed

reset!(m::model; profiles=true, img=false)

Erase existing profiles/raytrace status and invalidate the model's getVariable cache.

Call this after mutating ring fields directly (e.g. m.rings[1].I .= ...) so that subsequent getVariable/getProfile calls see the new values instead of memoized ones.

Parameters

• m::model: Model object to reset
• profiles::Bool=true: If true, reset profiles
• img::Bool=false: If true, reset raytracing boolean (does not change existing model but allows model to be raytraced again after combining other new models)

rotate3D(r::Float64, ϕ₀::Float64, i::Float64, rot::Float64, θₒ::Float64, reflect::Bool=false)

Transform from ring coordinates to 3D coordinates where camera is at +x.

Parameters

• r::Float64: Radius from central mass (in terms of rₛ)
• ϕ₀::Float64: Starting azimuthal angle in ring plane (rad)
• i::Float64: Inclination angle of ring plane (rad)
• rot::Float64: Rotation of system plane about z axis (rad)
• θₒ::Float64: Opening angle of point (rad)
• reflect::Bool=false: Whether to reflect across the ring plane

Returns

• Vector{Float64}: [x; y; z] coordinates in 3D space

Note

This method allocates its result; in performance-critical loops use rotate3D_scalar which returns a stack-allocated tuple instead.

rotate3D(r::Float64, ϕ₀::Float64, i::Float64, matrix::Matrix{Float64}, reflect::Bool=false)

Transform from ring coordinates to 3D coordinates where camera is at +x using a precomputed rotation matrix.

Parameters

• r::Float64: Radius from central mass (in terms of rₛ)
• ϕ₀::Float64: Starting azimuthal angle in ring plane (rad)
• i::Float64: Inclination angle of ring plane (rad)
• matrix::Matrix{Float64}: Precomputed rotation matrix for the ring
• reflect::Bool=false: Whether to reflect across the ring plane

Returns

• Vector{Float64}: [x; y; z] coordinates in 3D space

rotate3D_scalar(r::Float64, ϕ₀::Float64, i::Float64, rot::Float64, θₒ::Float64, reflect::Bool=false) -> NTuple{3,Float64}

Allocation-free version of rotate3D: transform from ring coordinates to 3D coordinates where camera is at +x, returning an (x, y, z) tuple instead of a heap-allocated Vector. Prefer this in hot loops.

rotate3D_vector_scalar(vx::Float64, vy::Float64, vz::Float64, i::Float64, rot::Float64, θₒ::Float64, reflect::Bool=false) -> NTuple{3,Float64}

Allocation-free rotation of an arbitrary vector (vx,vy,vz) from initial XY-plane coordinates into 3D space (camera at +x), with optional reflection across the ring plane. The matrix entries are those of get_r3D(i,rot,θₒ) inlined; used for rotating velocity vectors as well as positions.

profile(m::model, [variable], [kwargs...])
profile(rm::ResidentModel, variable, [kwargs...])
profile(cm::CompositeModel, [name], [kwargs...])
profile(rcm::ResidentCompositeModel, [name], [kwargs...])
profile(p::profile, [kwargs...])

Plot line profiles, normalized to the maximum value of each profile.

Arguments

• model: A model object containing profile data, or a ResidentModel.
• variable: A symbol or string (or list of symbols/strings) specifying which profile to plot. For a host model this is optional — if omitted, every profile already set on the model is plotted. A ResidentModel stores no preset profiles, so a variable is required there and the profile is computed on demand via getProfile (choose from :line, :delay, :r, :ϕ, :phase, :moment2).
• kwargs...: Additional keyword arguments passed to Plots.plot

CompositeModel and bare profile support (Workstream 4/5)

• cm::CompositeModel: plots the requested name profile (default :line, computed on demand via getProfile(cm, name; line=...)) for every line, one series per line labeled by the line name and normalized by that line's own max |binSums| (the same convention the multi-profile model case below uses) – so lines are shape-comparable regardless of fluxRatio. The one exception is profile(cm, :ratio, (a, b)) (W5): the velocity-resolved line ratio is computed from the PAIR of lines via getProfile(cm, :ratio; lines=(a, b)) and rendered as a single unnormalized series (the ratio's absolute scale is the point), identical to plotting the returned bare profile. The pair is a positional argument because Plots reserves lines as a magic line-attribute alias – a lines=... keyword would be consumed by the attribute preprocessing and never reach the recipe. Other :ratio keywords (bins, floor, ...) are not forwarded – for those, call getProfile yourself and plot the returned struct.
• rcm::ResidentCompositeModel: same as the CompositeModel case (one normalized series per line), except each line's profile is computed on its device-resident columns via getProfile(rcm, name; line=...) (W4-G3), with the resident restrictions (uniform bins, no custom Function profiles).
• p::profile: plots a single bare profile struct directly (binCenters vs binSums, unnormalized). Needed so Workstream 5's :ratio profile (also a plain profile struct) is directly plottable. Note profile is both this struct's type and (via @userplot Profile above) the name of the generated plotting function – they're the same generic function, but the @kwdef constructor takes zero positional arguments (kwargs only), while this plotting method takes one-or-more positional arguments, so profile(p) is never ambiguous with profile(name=..., binCenters=..., ...).

+(rm1::ResidentModel, rm2::ResidentModel) -> ResidentModel

Combine two device-resident models by concatenating their columns on their existing backend — when both live on the GPU the vcats stay on the device, so this is the fast path for gpu(m1) + gpu(m2) (or for merging submodels built on-device with residentDiskWindModel/residentCloudModel) without a host round-trip. Mirrors +(::model, ::model): submodels are stacked and nSubModels adds. Both handles must share the same backend type and element type.

Like the host +, this stacks submodels without raytracing, but it does preserve raytrace metadata when both operands carry it (merged via _rt_merge_meta). Use raytrace!(rm) for the device-resident τ-scan, or host raytrace!(m1 + m2; backend=…) for the one-shot host-orchestrated path.

Base.getindex(rcm::ResidentCompositeModel, line::String) -> ResidentModel

Retrieve the ResidentModel registered for line. Errors (listing the known line names) if line is not present, matching the host CompositeModel indexing behavior.

Base.length(rcm::ResidentCompositeModel) -> Int

Number of lines registered in rcm.

Base.show(io::IO, rcm::ResidentCompositeModel)

Print one row per line — name, lineCenter, fluxRatio, number of flattened points, element type, and backend — mirroring the host CompositeModel show.

flatten(m::model; T=Float64) -> ModelArrays{T,Vector{T},Vector{Bool}}

Build a flat structure-of-arrays snapshot of m, one entry per (ring, point), as the host-side bridge to the GPU (gpu(m) adapts the result to CuArray).

Each column is gathered through the package's existing, audited expansion machinery rather than a GPU-specific copy: per-point fields and per-ring scalars are aligned elementwise by expandPerPoint (so e.g. a scalar η is broadcast to every point of its ring), the 3D coordinates come from the ring getXYZ cache, and the camera columns from m.camera. The point ordering matches getVariable(m, sym, flatten=true).

getVariable(rm::ResidentModel, variable::Union{String,Symbol}) -> Vector

Host copy of a single flat ModelArrays column by name (e.g. :I, :v, :α). This is the recipe-facing counterpart of getVariable(::model, …); the columns are already flat so it takes no flatten keyword, and it does not accept a Function (compute custom quantities on the host model). When the columns live on the GPU the result is copied back to the host.

gpu(m::model; T=Float32) -> ResidentModel
gpu(ma::ModelArrays) -> ModelArrays

Move a host model or flat ModelArrays snapshot onto the GPU and return a device-resident handle for repeated observable calls. Requires CUDA.jl to be loaded so the package extension can provide the CUDA-backed methods.

resident(m::model; T=Float64, backend=KernelAbstractions.CPU(), raytrace=false) -> ResidentModel

Flatten m and wrap it as a ResidentModel on backend. The default CPU() backend keeps everything on the host — useful for testing the resident pipeline without a GPU. Use gpu(m) (with CUDA.jl loaded) to build a device-resident handle.

Pass raytrace=true to attach the metadata that raytrace!(::ResidentModel) needs. It is only built for host models with more than one submodel (single-submodel handles keep rt === nothing); on-device constructors carry their own metadata for later combination via +.

ModelArrays{T,V,B}

Flat structure-of-arrays representation of a model, with one entry per model point for geometry, velocity, intensity, response, camera coordinates, and reflection state. This is the common host/device payload used by resident, gpu, and the resident observable kernels.

ResidentCompositeModel

Host-side mirror of CompositeModel for device-resident multi-line work: the same lines/lineCenters/fluxRatios bookkeeping, but each line's model is a ResidentModel handle instead of a host model. Each line stays its own ResidentModel — there is no device-side composite kernel work, only orchestration (the resident _fluxWeights/getSpectrum/getProfile/ lineOverlap methods in gpu_observables.jl run per-line device reductions/kernels and combine the small results on the host; the image/plot3d forwarding and the profile/spectrum recipe support mirror the host composite's — see composite.jl W4-G3 and the recipe docstrings).

Built by resident(cm::CompositeModel) (CPU backend) or gpu(cm::CompositeModel) (device) — both of those constructors live in src/composite.jl, NOT here: this file is included before composite.jl (see src/BroadLineRegions.jl), so a method signature annotated ::CompositeModel would UndefVarError at load time. The struct itself only references ResidentModel and can live here with the rest of the resident machinery.

Fields

• lines::Vector{String}: unique line names, in creation order (lines[1] is the reference line, as on the host CompositeModel).
• models::Dict{String,ResidentModel}: the per-line device-resident handle.
• lineCenters::Dict{String,Float64}: central wavelength per line, any consistent units.
• fluxRatios::Dict{String,Float64}: each line's velocity-integrated flux relative to lines[1].

ResidentModel(ma::ModelArrays, backend, nSubModels::Int [, rt])

A model held resident on a compute device: a flat ModelArrays snapshot plus the KernelAbstractions backend its columns live on. Built once (gpu(m) on the device, or resident on the CPU) and reused across many observable calls — getProfile, getΨ, getΨt, phase, secondMoment have methods on it that run the GPU kernels without re-flattening or re-transferring the model each call.

The optional rt field carries device-resident raytrace metadata (output grid + per-point submodel / discrete info; see RaytraceMeta) so raytrace!(::ResidentModel) can run the bin→sort→scan→compact entirely on backend. It is nothing for handles that were not built for raytracing; the observable methods ignore it.

_build_cloud_modelarrays(nClouds, seed; μ, β, F, rₛ, θₒ, γ, ξ, i, rescale,
    ηₒ, η₁, αRM, rNorm, backend, T) -> ModelArrays

Allocate the device columns on backend and run _rt_build_cloud_kernel! to draw nClouds clouds (one thread each) into a ModelArrays.

_build_diskwind_modelarrays(rMin, rMax, inc, nr, nϕ, scale, rot, θₒ, Ifun, vfun;
    ηₒ, η₁, αRM, rNorm, backend, T) -> ModelArrays

Allocate the device columns on backend and run _rt_build_disk_kernel! to fill them, then wrap as a ModelArrays. Ifun/vfun are GPU-safe scalar callables (r, ϕ, inc) -> value.

_rt_disk_intensity_fn(::Type{T}, f1, f2, f3, f4, α) -> closure

GPU-safe scalar intensity closure for the built-in DiskWindIntensity, capturing the (T-typed) shape factors so the fused construction kernel can call Ifun(r, ϕ, inc). Mirrors DiskWind_I_(r, ϕ, i, f1, f2, f3, f4, α).

_rt_disk_radial_velocity_fn(::Type{T}, vᵣFrac, inflow, rₛ) -> closure

GPU-safe scalar velocity closure for the built-in vCircularRadialDisk (Keplerian + radial in/outflow). With vᵣFrac = 0 it is identical to vCircularDisk.

_rt_disk_velocity_fn(::Type{T}, rₛ) -> closure

GPU-safe scalar velocity closure for the built-in vCircularDisk, capturing the (T-typed) Schwarzschild radius so the fused construction kernel can call vfun(r, ϕ, inc).

gpuCloudModel(nClouds, seed; T=Float32, kwargs...) -> ResidentModel

Build a cloud model directly on the GPU (device-resident ResidentModel) — the CUDA counterpart of residentCloudModel, defaulting to Float32. Requires CUDA.jl loaded.

gpuDiskWindModel(rMin, rMax, i; T=Float32, kwargs...) -> ResidentModel
gpuDiskWindModel(r̄, rFac, α, i; T=Float32, kwargs...) -> ResidentModel

Build a DiskWind model directly on the GPU (device-resident ResidentModel) — the CUDA counterpart of residentDiskWindModel, defaulting to Float32 (GeForce FP64 is ~1/64 rate). Requires CUDA.jl loaded (defined by the CUDA package extension); errors otherwise.

residentCloudModel(nClouds, seed; μ=500.0, β=1.0, F=0.5, rₛ=1.0, θₒ=π/2, γ=1.0, ξ=1.0, i=0.0,
    rescale=1.0, ηₒ=0.5, η₁=0.5, αRM=0.0, rNorm=1.0,
    backend=KernelAbstractions.CPU(), T=Float64) -> ResidentModel

Build a Pancoast-style cloud model directly on backend as a device-resident ResidentModel — the on-device counterpart of cloudModel(nClouds; rng=:philox, seed=…). Each cloud is drawn from an independent counter-based (Philox4x32) substream keyed by (seed, cloud index), so the realization is deterministic, seed-reproducible, and identical across CPU/CUDA backends.

seed is required (the stream is counter-based).

Two built-in physics options are supported on-device (selected per call):

• intensity — isotropic (intensity=:isotropic, scaled by rescale) or the anisotropic cloudIntensity (intensity=:cloud, parameter κ).
• velocity — vCircularCloud (velocity=:circular) or the Pancoast vCloudTurbulentEllipticalFlow (velocity=:turbulent, parameters σρᵣ, σρc, σΘᵣ, σΘc, θₑ, fEllipse, fFlow, σₜ). The turbulent path draws extra normals/uniform per cloud from the same substream.

For physics beyond these built-ins, build on the host and call gpu(m).

residentDiskWindModel(r̄, rFac, α, i; rot=0.0, nr=128, nϕ=256, scale=:log, kwargs...) -> ResidentModel

r̄/rFac/α parameterization (matching DiskWindModel(r̄, rFac, α, i; ...)): derives rMin/rMax via get_rMinMaxDiskWind and sets the intensity power-law to α.

residentDiskWindModel(rMin, rMax, i; nr=128, nϕ=256, scale=:log, rot=0.0, θₒ=0.0,
    f1, f2, f3, f4, α, ηₒ=0.5, η₁=0.5, αRM=0.0, rNorm=1.0, rₛ=1.0,
    intensity=nothing, velocity=nothing, backend=KernelAbstractions.CPU(), T=Float64) -> ResidentModel

Build a DiskWind model directly on backend as a device-resident ResidentModel, without constructing host rings — the on-device counterpart of DiskWindModel(...) |> resident. One fused kernel generates the camera grid and computes geometry, velocity, intensity, and reverberation response per pixel. Reproduces flatten(DiskWindModel(rMin, rMax, i; ...)) to floating-point rounding (see the boundary-masking note in src/gpu_construct.jl).

The default physics is the built-in DiskWindIntensity (f1–f4, α) and vCircularDisk (rₛ). Pass vᵣFrac/inflow to switch the velocity to the built-in vCircularRadialDisk (Keplerian + radial in/outflow; vᵣFrac = 0 is vCircularDisk). Advanced users may instead pass intensity/velocity as GPU-safe scalar callables (r, ϕ, inc) -> value (e.g. a let-closure over isbits parameters) to run custom physics in the same fused kernel; generic custom Functions that cannot run in a kernel should use host construction + gpu(m) instead.

Use gpu-backed construction (CUDA loaded) via gpuDiskWindModel; the default CPU() backend builds on the host and is mainly for testing the on-device pipeline without a GPU.

_finiteVRange(rm::ResidentModel) -> (vmin, vmax)

Resident counterpart of _finiteVRange(::model): min/max of the device v column over points where isfinite(v) && isfinite(I*ΔA) (the W4-T4 finiteness rule — finite I*ΔA alone would let a NaN v poison the range, and vice versa), computed as two device reductions with ±Inf sentinels in the style of _rt_device_nanmin/_rt_device_nanmax. Returns host Float64s; (NaN, NaN) when the model has no finite points, exactly like the host helper.

_fluxWeights(rcm::ResidentCompositeModel) -> Dict{String,Float64}

Resident counterpart of _fluxWeights(::CompositeModel): per-line Σ I*ΔA as a single device reduction (result to host), counting only points where isfinite(v) && isfinite(I*ΔA) — the same finiteness rule as the CPU helper (both conditions, matching exactly what the histogram kernels count; a plain sum would return NaN the moment a sentinel survives). Returns fluxRatios[line] / totalFlux(line) per line.

Errors (naming the offending line) if any line's total flux is not positive, matching the CPU zero-flux guard: a zero-emission line cannot be normalized to unit integral, and an Inf weight would make the line silently vanish from the spectrum instead.

getProfile(rcm::ResidentCompositeModel, name; line=nothing, lines=nothing, kwargs...) -> profile

Per-line forwarding of the resident getProfile: compute name's profile for the ResidentModel registered under line (getProfile(rcm[line], name; kwargs...)), running the device kernels on that line's resident columns.

Mirrors getProfile(::CompositeModel, …) exactly – a single method covering the whole positional signature, branching at runtime (Julia cannot dispatch on keywords):

• name === :ratio (velocity-resolved line ratio, W5-G1) requires lines::Tuple{String,String} and computes the ratio of the two lines' device-binned velocity profiles with identical keyword semantics to the host method (see the :ratio section of getProfile(::CompositeModel, …) for the full argument docs): shared bins from the union of the two lines' _finiteVRange device reductions, one device binned sum per line, the division and NaN/floor handling on the small host vectors via the shared _ratioDivide. The resident restriction applies: a user-supplied bins edge vector must be uniform (_rt_resident_edges).
• every other name requires the line keyword and forwards to the per-line ResidentModel method, whose restrictions apply (uniform bins only, no custom dx, no custom-Function profiles).

getProfile(rm::ResidentModel, name; bins=100, dx=nothing, kwargs...) -> profile

GPU-resident getProfile. Supports :line, :delay, :r, :ϕ, :phase, and :moment2, matching the host getProfile(::model, …) semantics (binnedSum's left-exclusive bins). Custom-Function profiles and a custom dx integration element are not supported on the resident path (the model's ΔA is always used); pass a host model for either.

getSpectrum(rcm::ResidentCompositeModel; bins=100, z=0.0, minX=nothing, maxX=nothing,
            centered=false, overflow=nothing, kwargs...)
    -> (edges, centers, flux::Dict{String,Vector{Float64}}, total::Vector{Float64})

Resident counterpart of getSpectrum(::CompositeModel), with an identical return shape and identical keyword semantics (see the CPU method for the full argument docs). The shared wavelength range is computed on the host from per-line device min/max reductions of v (_finiteVRange, NaN-aware) unless pinned by the caller via minX/maxX; each line is then binned on its own device via the existing resident binned-sum path with the shared bins/minX/maxX (constructBinEdges is deterministic, so every line gets bit-identical edges), and the per-line host vectors are combined into total.

As on the CPU, the binning keywords (minX/maxX/centered/overflow) forward to every line's binned sum, and overflow defaults to true for integer bins (flux-preserving – non-centered edges built at exactly [λmin, λmax] place each line's extremal points on the boundary, which the binning otherwise drops; a user-narrowed window accumulates out-of-window flux into the boundary bins) and false for a user-supplied edge vector; an explicit overflow always wins. The resident path additionally requires a user edge vector to be uniform (the device histogram kernels assume uniform spacing; see _rt_resident_edges).

getΨ(rm::ResidentModel, vEdges::AbstractVector, tEdges::AbstractVector) -> Ψ

GPU-resident 2D transfer function Ψ(v, t); uniform edges only. Matches the host getΨ.

getΨt(rm::ResidentModel, tEdges::AbstractVector, overflow::Bool=false) -> Ψt

GPU-resident 1D transfer function Ψ(t); uniform tEdges only. Matches the host getΨt (binnedSum's left-exclusive interior-edge convention, empty/NaN-poisoned bins floored to 1e-30, optional overflow into the edge bins).

lineOverlap(rcm::ResidentCompositeModel) -> Vector{NamedTuple}

Resident counterpart of lineOverlap with identical pairwise-overlap semantics, via the shared _lineOverlap implementation (src/composite.jl); the only difference is that each line's velocity range comes from the device reductions of _finiteVRange(::ResidentModel) (W4-G2) instead of the host gathers. This is also what lets the spectrum plot recipe accept a ResidentCompositeModel with no recipe changes.

lineRatio(rcm::ResidentCompositeModel, a::String, b::String) -> Float64

Resident counterpart of lineRatio via the shared _lineRatio implementation (src/composite.jl) – fluxRatios is host metadata on both composite types, so no device work is involved.

phase(rm::ResidentModel; U, V, PA, BLRAng, returnAvg=false, offAxisInds=nothing, bins=100, kwargs...)

GPU-resident differential phase, matching phase(::model, …).

secondMoment(rm::ResidentModel; U, V, PA, BLRAng, returnAvg=false, offAxisInds=nothing, bins=100, kwargs...)

GPU-resident image second moment, matching secondMoment(::model, …).